From beaf4c59137554942259a40d9c1ff3aa26bd76c4 Mon Sep 17 00:00:00 2001 From: Will Simmons Date: Fri, 8 Nov 2024 14:42:47 +0000 Subject: [PATCH 1/6] Copy files over to a new branch just for this feature --- tket/CMakeLists.txt | 3 + tket/include/tket/Clifford/APState.hpp | 146 ++++ tket/include/tket/Converters/Converters.hpp | 31 + tket/src/Clifford/APState.cpp | 861 ++++++++++++++++++++ tket/src/Converters/APStateConverters.cpp | 176 ++++ tket/test/CMakeLists.txt | 1 + tket/test/src/test_APState.cpp | 347 ++++++++ 7 files changed, 1565 insertions(+) create mode 100644 tket/include/tket/Clifford/APState.hpp create mode 100644 tket/src/Clifford/APState.cpp create mode 100644 tket/src/Converters/APStateConverters.cpp create mode 100644 tket/test/src/test_APState.cpp diff --git a/tket/CMakeLists.txt b/tket/CMakeLists.txt index 41808ece11..e7cfc3a681 100644 --- a/tket/CMakeLists.txt +++ b/tket/CMakeLists.txt @@ -168,9 +168,11 @@ target_sources(tket src/Circuit/SubcircuitFinder.cpp src/Circuit/ThreeQubitConversion.cpp src/Circuit/ToffoliBox.cpp + src/Clifford/APState.cpp src/Clifford/ChoiMixTableau.cpp src/Clifford/SymplecticTableau.cpp src/Clifford/UnitaryTableau.cpp + src/Converters/APStateConverters.cpp src/Converters/ChoiMixTableauConverters.cpp src/Converters/Gauss.cpp src/Converters/PauliGraphConverters.cpp @@ -322,6 +324,7 @@ target_sources(tket include/tket/Circuit/StatePreparation.hpp include/tket/Circuit/ThreeQubitConversion.hpp include/tket/Circuit/ToffoliBox.hpp + include/tket/Clifford/APState.hpp include/tket/Clifford/ChoiMixTableau.hpp include/tket/Clifford/SymplecticTableau.hpp include/tket/Clifford/UnitaryTableau.hpp diff --git a/tket/include/tket/Clifford/APState.hpp b/tket/include/tket/Clifford/APState.hpp new file mode 100644 index 0000000000..402106b11c --- /dev/null +++ b/tket/include/tket/Clifford/APState.hpp @@ -0,0 +1,146 @@ +// Copyright 2019-2024 Cambridge Quantum Computing +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#pragma once + +#include "tket/OpType/OpType.hpp" +#include "tket/Utils/MatrixAnalysis.hpp" + +namespace tket { + +/** + * APState class for Clifford states with global phase tracking + * + * The "affine with phases" form of a Clifford state from ZX calculus (see + * Kissinger & van de Wetering, "Picturing Quantum Software") represents n-qubit + * Clifford states uniquely with the following data: + * - A binary (k,n) matrix A in reduced row echelon form. The leading columns of + * each row are referred to as "leading qubits", and the others as "free + * qubits". We apply H gates to all free qubits, then for each row of A we apply + * a CX targeted on the leading qubit and control by each other entry in the + * row. + * - A binary (k) vector B describing X gates applied to the leading qubits. + * - A Clifford phase polynomial Phi over the free qubits. Wlog this can be in + * the form of a symmetric, zero-diagonal, binary (n-k,n-k) matrix E describing + * CZ gates and a (n-k) vector P of integers mod 4 describing S gates. + * + * This gives a canonical statevector: + * \sum_{x, Ax=B} i^{\Phi(x)} \ket{x} + * + * This canonical statevector fixes a reference phase from which we can track + * global phase with an additional parameter. + * + * When applying gates, qubits may switch between being leading and free + * (including those that aren't involved in the gate due to the need to reduce + * to the canonical form, e.g. reduced row echelon form for A). The updates are + * easiest to prove in the ZX calculus for the gateset (CZ, S i.e. pi/2 green + * phase, SX i.e. pi/2 red phase). + * + * To save on resizing the matrices and vectors, we will keep them at full size + * and just assert correctness that, e.g. every entry in A is either on the + * diagonal (to indicate the qubit is leading) or in the row of a leading qubit + * and a column of a later following qubit. + */ +class APState { + public: + /** + * An upper triangular matrix whose diagonal entries indicate leading qubits + * and off-diagonal entries represent a CX gate from the column qubit + * (necessarily a free qubit) to the row qubit. + * + * If a diagonal entry is 0 (it is a free qubit), then every entry in the row + * is 0. If a diagonal entry is 1 (it is a leading qubit), all other entries + * in the column are 0. + */ + MatrixXb A; + + /** + * A vector indicating X gates on leading qubits. + * + * Must be 0 on every free qubit. + */ + VectorXb B; + + /** + * A symmetric, zero diagonal matrix whose entries indicate CZs between free + * qubits. + * + * Every diagonal entry must be 0. + * The column and row for each leading qubit must be 0. + */ + MatrixXb E; + + /** + * A vector indicating S^{P(i)} on qubit i which must be free. + * + * Must be 0 on every leading qubit. + */ + Eigen::VectorXi P; + + /** + * The global phase term (in half-turns). + */ + Expr phase; + + /** + * Constructs a state in AP form from given data. + */ + APState( + const MatrixXb& A_, const VectorXb& B_, const MatrixXb& E_, + const Eigen::VectorXi& P_, const Expr& phase_); + + /** + * Constructs the state \ket{0}^{\otimes n_qubits} in AP form. + */ + APState(unsigned n_qubits); + + /** + * Constructs the state in AP form from a given statevector. + */ + APState(const Eigen::VectorXcd& sv); + + /** + * Verifies the internal correctness of the data structure. Throws an + * exception if the data structure is invalid. + */ + void verify(); + + /** + * Calculates the statevector of the state. + */ + Eigen::VectorXcd to_statevector(); + + /** + * Applies a CZ gate to the state. O(n^2) wrt number of qubits in the state. + */ + void apply_CZ(unsigned ctrl, unsigned trgt); + + /** + * Applies an S gate to the state. O(n^2) wrt number of qubits in the state. + */ + void apply_S(unsigned q); + + /** + * Applies a V gate to the state. O(n^2) wrt number of qubits in the state. + */ + void apply_V(unsigned q); + + /** + * Applies an unparameterised Clifford gate to the state on the chosen qubits. + * O(n^2) wrt number of qubits in the state. + */ + void apply_gate(OpType type, const std::vector& qbs); +}; + +} // namespace tket \ No newline at end of file diff --git a/tket/include/tket/Converters/Converters.hpp b/tket/include/tket/Converters/Converters.hpp index a3de8c9146..0a9f9809cc 100644 --- a/tket/include/tket/Converters/Converters.hpp +++ b/tket/include/tket/Converters/Converters.hpp @@ -15,6 +15,7 @@ #pragma once #include "tket/Circuit/Circuit.hpp" +#include "tket/Clifford/APState.hpp" #include "tket/Clifford/ChoiMixTableau.hpp" #include "tket/Clifford/UnitaryTableau.hpp" #include "tket/PauliGraph/PauliGraph.hpp" @@ -45,6 +46,20 @@ Circuit unitary_rev_tableau_to_circuit(const UnitaryRevTableau &tab); */ ChoiMixTableau circuit_to_cm_tableau(const Circuit &circ); +/** + * Construct an APState for a given circuit. + * Assumes all input qubits are initialised in the |0> state. + * Will throw an exception if it contains non-Clifford gates. + */ +APState circuit_to_apstate(const Circuit &circ); + +/** + * Construct a circuit producing the state described by the APState. + * Uses the standard circuit form implied by the ZX description of reduced + * AP-form (layered as X-H-CX-CZ-S). + */ +Circuit apstate_to_circuit(const APState &ap); + /** * Constructs a circuit producing the same effect as a ChoiMixTableau. * Since Circuit does not support distinct qubit addresses for inputs and @@ -130,6 +145,22 @@ std::pair cm_tableau_to_unitary_extension_circuit( const std::vector &post_names = {}, CXConfigType cx_config = CXConfigType::Snake); +/** + * Convert a SymplecticTableau describing a stabiliser state to reduced AP-form + * (i.e. an APState object). Since tableau methods do not track global phase, we + * set it to 0. Throws an exception if the tableau does not describe a pure + * state (i.e. it must have n commuting rows for n qubits). + */ +APState tableau_to_apstate(SymplecticTableau tab); + +/** + * Convert an APState describing a stabiliser state into a tableau form, + * discarding the global phase information. This allows us to reuse tableau + * synthesis methods to synthesise APState objects (the result would need + * re-simulating as an APState to compare and deduce the required global phase). + */ +SymplecticTableau apstate_to_tableau(const APState &ap); + PauliGraph circuit_to_pauli_graph(const Circuit &circ); /** diff --git a/tket/src/Clifford/APState.cpp b/tket/src/Clifford/APState.cpp new file mode 100644 index 0000000000..a21e4b446b --- /dev/null +++ b/tket/src/Clifford/APState.cpp @@ -0,0 +1,861 @@ +// Copyright 2019-2024 Cambridge Quantum Computing +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include "tket/Clifford/APState.hpp" + +#include + +#include "tket/OpType/OpTypeInfo.hpp" + +namespace tket { + +APState::APState( + const MatrixXb& A_, const VectorXb& B_, const MatrixXb& E_, + const Eigen::VectorXi& P_, const Expr& phase_) + : A(A_), B(B_), E(E_), P(P_), phase(phase_) { + verify(); +} + +APState::APState(unsigned n_qubits) + : A(MatrixXb::Identity(n_qubits, n_qubits)), + B(VectorXb::Zero(n_qubits)), + E(MatrixXb::Zero(n_qubits, n_qubits)), + P(Eigen::VectorXi(n_qubits)), + phase(0) {} + +unsigned clifford_phase(const Complex& c) { + unsigned res = 0; + unsigned n_possibles = 0; + if (c.real() > EPS) { + res = 0; + ++n_possibles; + } + if (c.imag() > EPS) { + res = 1; + ++n_possibles; + } + if (c.real() < -EPS) { + res = 2; + ++n_possibles; + } + if (c.imag() < -EPS) { + res = 3; + ++n_possibles; + } + TKET_ASSERT(n_possibles == 1); + return res; +} + +APState::APState(const Eigen::VectorXcd& sv) { + unsigned n_qbs = 0; + while (sv.size() > (1 << n_qbs)) ++n_qbs; + TKET_ASSERT(sv.size() == (1 << n_qbs)); + + // Find non-zero entries as a vector space and offset + unsigned z0 = 0; + for (unsigned x = 0; x < sv.size(); ++x) { + if (sv(x) != 0.) { + z0 = x; + break; + } + } + std::vector offsets; + unsigned n_non_zero = 0; + for (unsigned x = 1; x < sv.size(); ++x) { + if (sv(z0 ^ x) != 0.) { + ++n_non_zero; + if (n_non_zero == (1 << offsets.size())) offsets.push_back(x); + } + } + + // Find A as the dual space + MatrixXb offset_mat(offsets.size(), n_qbs); + for (unsigned r = 0; r < offsets.size(); ++r) { + unsigned off = offsets.at(r); + for (unsigned c = 0; c < n_qbs; ++c) { + // Binary encoding of offsets in reverse order to guarantee free qubits + // are the later ones, meaning we produce A in row echelon form + offset_mat(r, c) = (off & (1 << c)) != 0; + } + } + std::vector> row_ops = + gaussian_elimination_row_ops(offset_mat); + for (const std::pair& op : row_ops) { + for (unsigned j = 0; j < n_qbs; ++j) { + offset_mat(op.second, j) ^= offset_mat(op.first, j); + } + } + std::map mat_leaders; + A = MatrixXb::Zero(n_qbs, n_qbs); + for (unsigned c = 0; c < n_qbs; ++c) { + bool free_qubit = false; + for (unsigned r = 0; r < offsets.size(); ++r) { + if (offset_mat(r, c)) { + auto inserted = mat_leaders.insert({r, c}); + if (inserted.second) { + free_qubit = true; + break; + } else { + // Reverse bit orderings back to normal + A(n_qbs - 1 - c, n_qbs - 1 - inserted.first->second) = true; + } + } + } + A(n_qbs - 1 - c, n_qbs - 1 - c) = !free_qubit; + } + + VectorXb z0_vec(n_qbs); + for (unsigned i = 0; i < n_qbs; ++i) { + z0_vec(i) = ((z0 & (1 << (n_qbs - 1 - i))) != 0); + } + B = A * z0_vec; + + E = MatrixXb::Zero(n_qbs, n_qbs); + P = Eigen::VectorXi::Zero(n_qbs); + unsigned neutral_z0 = z0; // Index with 0 for all free qubits + std::map offset_for_free; + for (const std::pair& row_qfree : mat_leaders) { + unsigned offset = 0; + for (unsigned i = 0; i < n_qbs; ++i) + if (offset_mat(row_qfree.first, i)) offset += (1 << i); + unsigned qfree = n_qbs - 1 - row_qfree.second; + offset_for_free.insert({qfree, offset}); + if ((neutral_z0 & (1 << row_qfree.second)) != 0) + neutral_z0 = neutral_z0 ^ offset; + } + for (const std::pair& qfree_offset : + offset_for_free) { + unsigned qfree = qfree_offset.first; + unsigned offset = qfree_offset.second; + auto complex_phase = sv(neutral_z0 ^ offset) / sv(neutral_z0); + P(qfree) = clifford_phase(complex_phase); + for (const std::pair& qfree_offset2 : + offset_for_free) { + if (qfree_offset == qfree_offset2) break; // Only go up to solved phases + unsigned qfree2 = qfree_offset2.first; + unsigned offset2 = qfree_offset2.second; + auto complex_phase_cz = + sv(neutral_z0 ^ offset ^ offset2) / sv(neutral_z0); + E(qfree, qfree2) = E(qfree2, qfree) = + ((clifford_phase(complex_phase_cz) - P(qfree) - P(qfree2)) % 4) == 2; + } + } + + phase = Expr(std::arg(sv(neutral_z0)) / PI); +} + +void APState::verify() { + // Check dimensions all agree + unsigned n_qubits = A.rows(); + TKET_ASSERT(A.cols() == n_qubits); + TKET_ASSERT(B.size() == n_qubits); + TKET_ASSERT(E.rows() == n_qubits); + TKET_ASSERT(E.cols() == n_qubits); + TKET_ASSERT(P.size() == n_qubits); + + for (unsigned r = 0; r < n_qubits; ++r) { + for (unsigned c = 0; c < r; ++c) { + // Check A is upper triangular + TKET_ASSERT(!A(r, c)); + // Check E is symmetric + TKET_ASSERT(E(r, c) == E(c, r)); + } + if (A(r, r)) { + // Check leaders are eliminated in A (row echelon) + for (unsigned r2 = 0; r2 < r; ++r2) { + TKET_ASSERT(!A(r2, r)); + } + } else { + // Check A only has entries in rows with a leader + for (unsigned c = r + 1; c < n_qubits; ++c) { + TKET_ASSERT(!A(r, c)); + } + } + } + + for (unsigned q = 0; q < n_qubits; ++q) { + if (A(q, q)) { + // Check P is zero on leaders + TKET_ASSERT(P(q) % 4 == 0); + // Check E is zero on leaders (simplified by symmetry) + for (unsigned q2 = 0; q2 < n_qubits; ++q2) { + TKET_ASSERT(!E(q, q2)); + } + } else { + // Check B is zero on free qubits + TKET_ASSERT(!B(q)); + // Check diagonals of E are still zero + TKET_ASSERT(!E(q, q)); + } + } +} + +VectorXb z2_mult(const MatrixXb& M, const VectorXb& v) { + VectorXb res = VectorXb::Zero(M.rows()); + for (unsigned i = 0; i < M.cols(); ++i) { + if (v(i)) { + for (unsigned j = 0; j < M.rows(); ++j) res(j) ^= M(j, i); + } + } + return res; +} + +Eigen::VectorXcd APState::to_statevector() { + unsigned n_qubits = A.cols(); + Eigen::VectorXcd sv = Eigen::VectorXcd::Zero(1 << n_qubits); + unsigned n_terms = 0; + Complex g_phase = std::exp(i_ * PI * eval_expr(phase).value()); + for (unsigned x = 0; x < 1 << n_qubits; ++x) { + VectorXb x_binary = VectorXb::Zero(n_qubits); + for (unsigned i = 0; i < n_qubits; ++i) { + unsigned mask = 1 << (n_qubits - 1 - i); + x_binary(i) = ((x & mask) != 0); + } + + if (z2_mult(A, x_binary) == B) { + ++n_terms; + unsigned i_phases = 0; + for (unsigned q = 0; q < n_qubits; ++q) { + if (x_binary(q)) { + i_phases += P(q); + for (unsigned q2 = q; q2 < n_qubits; ++q2) { + if (E(q, q2) && x_binary(q2)) i_phases += 2; + } + } + } + switch (i_phases % 4) { + case 0: { + sv(x) = g_phase; + break; + } + case 1: { + sv(x) = i_ * g_phase; + break; + } + case 2: { + sv(x) = -g_phase; + break; + } + case 3: { + sv(x) = -i_ * g_phase; + break; + } + default: { + TKET_ASSERT(false); + } + } + } + } + return pow(n_terms, -0.5) * sv; +} + +void APState::apply_CZ(unsigned ctrl, unsigned trgt) { + if (ctrl == trgt) + throw std::logic_error( + "APState::apply_CZ: cannot apply CZ when control and target qubits " + "match"); + unsigned n_qubits = A.cols(); + if (ctrl >= n_qubits || trgt >= n_qubits) + throw std::logic_error("APState::apply_CZ: qubit indices out of range"); + + if (A(ctrl, ctrl)) { + if (A(trgt, trgt)) { + // ctrl and trgt are both leading + + std::set just_ctrl, just_trgt, both; + // Add phases and determine connectivity + for (unsigned q = 0; q < n_qubits; ++q) { + if (A(ctrl, q)) { + if (B(trgt)) P(q) += 2; + if (A(trgt, q)) { + if (!B(ctrl)) P(q) += 2; + both.insert(q); + } else { + just_ctrl.insert(q); + } + } else if (A(trgt, q)) { + if (B(ctrl)) P(q) += 2; + just_trgt.insert(q); + } + } + // ctrl and trgt were included in the previous loop, so reset them + just_ctrl.erase(ctrl); + just_trgt.erase(trgt); + P(ctrl) = 0; + P(trgt) = 0; + // Pivoting-style update between those connected to ctrl, trgt, and both + for (const unsigned q1 : just_ctrl) { + for (const unsigned q2 : just_trgt) { + E(q1, q2) ^= true; + E(q2, q1) ^= true; + } + for (const unsigned q2 : both) { + E(q1, q2) ^= true; + E(q2, q1) ^= true; + } + } + for (const unsigned q1 : just_trgt) { + for (const unsigned q2 : both) { + E(q1, q2) ^= true; + E(q2, q1) ^= true; + } + } + + if (B(ctrl) && B(trgt)) phase += 1; + } else { + // ctrl is leading, trgt is free + + // Flip connectivity between trgt and every q on ctrl + for (unsigned q = ctrl + 1; q < n_qubits; ++q) { + if (A(ctrl, q)) { + E(trgt, q) ^= true; + E(q, trgt) ^= true; + // Note that if q = trgt, this does nothing, preserving the zero + // diagonal of E + } + } + // Add phase to trgt + if (A(ctrl, trgt) ^ B(ctrl)) P(trgt) += 2; + + // No global phase change + } + } else { + if (A(trgt, trgt)) { + // trgt is leading, ctrl is free + + // CZ is symmetric, so reuse the previous case + apply_CZ(trgt, ctrl); + } else { + // ctrl and trgt are both free + + // Just add the CZ gate to E + E(ctrl, trgt) ^= true; + E(trgt, ctrl) ^= true; + + // No global phase change + } + } +} + +void APState::apply_S(unsigned q) { + unsigned n_qubits = A.cols(); + if (q >= n_qubits) + throw std::logic_error("APState::apply_S: qubit index out of range"); + + if (A(q, q)) { + // q is leading + + // Local complementation update to E + for (unsigned q1 = q + 1; q1 < n_qubits; ++q1) { + if (A(q, q1)) { + for (unsigned q2 = q + 1; q2 < n_qubits; ++q2) { + E(q1, q2) ^= ((q1 != q2) && A(q, q2)); + } + } + } + + // Update global phase + if (B(q)) phase += .5; + + // Update local phases + unsigned local_phase_change = B(q) ? 3 : 1; + for (unsigned q1 = q + 1; q1 < n_qubits; ++q1) { + if (A(q, q1)) P(q1) += local_phase_change; + } + } else { + // q is free + + // Add to local phase + P(q) += 1; + + // No global phase change + } +} + +void APState::apply_V(unsigned q) { + unsigned n_qubits = A.cols(); + if (q >= n_qubits) + throw std::logic_error("APState::apply_V: qubit index out of range"); + + if (A(q, q)) { + // q is leading, but will become free + + // Local complementation for E + for (unsigned q1 = q; q1 < n_qubits; ++q1) { + if (A(q, q1)) { + for (unsigned q2 = q1 + 1; q2 < n_qubits; ++q2) { + E(q1, q2) ^= A(q, q2); + E(q2, q1) ^= A(q, q2); + } + + // Local phases + P(q1) += B(q) ? 1 : 3; + } + } + + // Global phase update + if (B(q)) phase += -.5; + + // q is no longer a leader + for (unsigned q1 = q; q1 < n_qubits; ++q1) { + A(q, q1) = false; + } + B(q) = false; + } else { + // q is free + + std::optional last_connected_leader; + for (unsigned q1 = q; q1 > 0;) { + --q1; + if (A(q1, q)) { + last_connected_leader = q1; + break; + } + } + + if (last_connected_leader.has_value()) { + // q is free and connected to a leader + + // For each other leader connected to q, add LCL's row in A + std::list other_leaders; + for (unsigned q1 = 0; q1 < *last_connected_leader; ++q1) { + if (A(q1, q)) { + other_leaders.push_back(q1); + for (unsigned q2 = *last_connected_leader; q2 < n_qubits; ++q2) { + A(q1, q2) ^= A(*last_connected_leader, q2); + } + } + } + + // Sort the connected frees into just q, just LCL, and both + std::list just_q, just_lcl, both; + for (unsigned q1 = 0; q1 < n_qubits; ++q1) { + if (q1 == q || q1 == *last_connected_leader) continue; + if (E(q, q1)) { + if (A(*last_connected_leader, q1)) + both.push_back(q1); + else + just_q.push_back(q1); + } else if (A(*last_connected_leader, q1)) + just_lcl.push_back(q1); + } + + // Split case by phase on q + // In each case: + // - Perform appropriate complementations between connected frees + // - Set LCL's connectivity in E appropriately + // - Modify q's connectivity in E + // - Add all local and global phases by a case matrix + if (P(q) % 2 == 0) { + // Handle complementations between neighbours + for (const unsigned q1 : just_lcl) { + // Local complement within group + for (const unsigned q2 : just_lcl) { + E(q1, q2) ^= true; + } + E(q1, q1) ^= true; + // Complement between groups + for (const unsigned q2 : just_q) { + E(q1, q2) ^= true; + E(q2, q1) ^= true; + } + // Add connections with q + E(q, q1) ^= true; + E(q1, q) ^= true; + } + for (const unsigned q1 : both) { + // Local complement within group + for (const unsigned q2 : both) { + E(q1, q2) ^= true; + } + E(q1, q1) ^= true; + // Complement between groups + for (const unsigned q2 : just_q) { + E(q1, q2) ^= true; + E(q2, q1) ^= true; + } + } + // Remove connections with q + for (const unsigned q1 : just_q) { + E(q, q1) ^= true; + E(q1, q) ^= true; + } + + // Move LCL from A to E + E(q, *last_connected_leader) = true; + E(*last_connected_leader, q) = true; + for (const unsigned q1 : just_q) { + E(q1, *last_connected_leader) = true; + E(*last_connected_leader, q1) = true; + } + for (const unsigned q1 : just_lcl) { + E(q1, *last_connected_leader) = true; + E(*last_connected_leader, q1) = true; + A(*last_connected_leader, q1) = false; + } + for (const unsigned q1 : both) A(*last_connected_leader, q1) = false; + A(*last_connected_leader, q) = false; + A(*last_connected_leader, *last_connected_leader) = false; + + // Phases + bool a = (P(q) % 4 == 2); + bool b = B(*last_connected_leader); + P(q) = b ? 1 : 3; + P(*last_connected_leader) = (a ^ b) ? 1 : 3; + B(*last_connected_leader) = false; + for (const unsigned q1 : other_leaders) B(q1) ^= b; + for (const unsigned q1 : both) P(q1) += a ? 3 : 1; + for (const unsigned q1 : just_lcl) P(q1) += (b ^ a) ? 1 : 3; + for (const unsigned q1 : just_q) P(q1) += b ? 2 : 0; + if (b) phase += a ? .5 : -.5; + } else { + // Handle complementations between neighbours + for (const unsigned q1 : just_lcl) { + // Complement between groups + for (const unsigned q2 : just_q) { + E(q1, q2) ^= true; + E(q2, q1) ^= true; + } + for (const unsigned q2 : both) { + E(q1, q2) ^= true; + E(q2, q1) ^= true; + } + // Add connections with q + E(q, q1) = true; + E(q1, q) = true; + } + for (const unsigned q1 : just_q) { + // Complement between groups + for (const unsigned q2 : both) { + E(q1, q2) ^= true; + E(q2, q1) ^= true; + } + // Remove connections with q + E(q, q1) = false; + E(q1, q) = false; + } + + // Move LCL from A to E + E(q, *last_connected_leader) = true; + E(*last_connected_leader, q) = true; + for (const unsigned q1 : just_q) { + E(q1, *last_connected_leader) = true; + E(*last_connected_leader, q1) = true; + } + for (const unsigned q1 : both) { + E(q1, *last_connected_leader) = true; + E(*last_connected_leader, q1) = true; + A(*last_connected_leader, q1) = false; + } + for (const unsigned q1 : just_lcl) + A(*last_connected_leader, q1) = false; + A(*last_connected_leader, q) = false; + A(*last_connected_leader, *last_connected_leader) = false; + + // Phases + bool a = (P(q) % 4 == 3); + bool b = B(*last_connected_leader); + P(q) = b ? 1 : 3; + P(*last_connected_leader) = a ? 2 : 0; + B(*last_connected_leader) = false; + for (const unsigned q1 : other_leaders) B(q1) ^= b; + for (const unsigned q1 : just_q) P(q1) += b ? 2 : 0; + for (const unsigned q1 : just_lcl) P(q1) += a ? 2 : 0; + for (const unsigned q1 : both) P(q1) += (a ^ b) ? 0 : 2; + if (a && b) phase += 1.; + } + } else if (P(q) % 2 == 0) { + // q has phase 0/pi and no connected leader + + // Local complementation update to E + std::list connected; + for (unsigned q1 = 0; q1 < n_qubits; ++q1) { + if (E(q, q1)) { + for (const unsigned q2 : connected) { + E(q1, q2) ^= true; + E(q2, q1) ^= true; + } + connected.push_back(q1); + + // Add local phase + P(q1) += (P(q) % 4 == 0) ? 1 : 3; + } + } + + // Local phase on q remains the same + + // Global phase + phase += (P(q) % 4 == 0) ? -0.25 : 0.25; + } else { + // q has phase +-pi/2 and no connected leader + + std::optional first_connected_free = std::nullopt; + for (unsigned q1 = 0; q1 < q; ++q1) { + if (E(q, q1)) { + first_connected_free = q1; + break; + } + } + + if (first_connected_free.has_value()) { + // q has phase +-pi/2, no connected leader, but a previous free which + // will become the new leader + + // Make FCF a leader + for (unsigned q1 = *first_connected_free; q1 < n_qubits; ++q1) { + if (E(q, q1)) { + A(*first_connected_free, q1) = true; + // q ends up disconnected in E + E(q, q1) = false; + E(q1, q) = false; + } + } + A(*first_connected_free, q) = true; + B(*first_connected_free) = (P(q) % 4 == 3); + + // Set phase of q to -pi/2 + P(q) = 3; + + // No global phase change + + // Make A row reduced + for (unsigned q1 = 0; q1 < *first_connected_free; ++q1) { + if (A(q1, *first_connected_free)) { + for (unsigned q2 = *first_connected_free; q2 < n_qubits; ++q2) { + A(q1, q2) ^= A(*first_connected_free, q2); + } + B(q1) ^= B(*first_connected_free); + } + } + + // Need to apply rules for S and CZ gates on FCF to remove + unsigned s_gates = P(*first_connected_free); + P(*first_connected_free) = 0; + std::list cz_targets; + E(*first_connected_free, q) = false; + E(q, *first_connected_free) = false; + for (unsigned q1 = 0; q1 < n_qubits; ++q1) { + if (E(*first_connected_free, q1)) { + cz_targets.push_back(q1); + E(*first_connected_free, q1) = false; + E(q1, *first_connected_free) = false; + } + } + for (unsigned i = 0; i < s_gates; ++i) apply_S(*first_connected_free); + for (const unsigned trgt : cz_targets) + apply_CZ(*first_connected_free, trgt); + } else { + // q has phase +-pi/2, no connected leader, and no previous free, so q + // will become the new leader + A(q, q) = true; + + std::list connected; + for (unsigned q1 = q + 1; q1 < n_qubits; ++q1) { + if (E(q, q1)) { + // Local complementation + for (const unsigned q2 : connected) { + E(q1, q2) ^= true; + E(q2, q1) ^= true; + } + connected.push_back(q1); + + // Add to A + A(q, q1) = true; + + // Reset E rows + E(q, q1) = false; + E(q1, q) = false; + } + } + + // Update local and global phases + if (P(q) % 4 == 1) { + for (const unsigned q1 : connected) P(q1) += 3; + } else { + for (const unsigned q1 : connected) P(q1) += 1; + B(q) = true; + phase += -.5; + } + P(q) = 0; + } + } + } +} + +void APState::apply_gate(OpType type, const std::vector& qbs) { + switch (type) { + case OpType::Z: { + apply_S(qbs.at(0)); + apply_S(qbs.at(0)); + break; + } + case OpType::X: { + apply_V(qbs.at(0)); + apply_V(qbs.at(0)); + phase += .5; + break; + } + case OpType::Y: { + apply_S(qbs.at(0)); + apply_S(qbs.at(0)); + apply_V(qbs.at(0)); + apply_V(qbs.at(0)); + phase += 1.; + break; + } + case OpType::S: { + apply_S(qbs.at(0)); + break; + } + case OpType::Sdg: { + apply_S(qbs.at(0)); + apply_S(qbs.at(0)); + apply_S(qbs.at(0)); + break; + } + case OpType::V: { + apply_V(qbs.at(0)); + break; + } + case OpType::SX: { + apply_V(qbs.at(0)); + phase += .25; + break; + } + case OpType::Vdg: { + apply_V(qbs.at(0)); + apply_V(qbs.at(0)); + apply_V(qbs.at(0)); + phase += 1.; + break; + } + case OpType::SXdg: { + apply_V(qbs.at(0)); + apply_V(qbs.at(0)); + apply_V(qbs.at(0)); + phase += .75; + break; + } + case OpType::H: { + apply_S(qbs.at(0)); + apply_V(qbs.at(0)); + apply_S(qbs.at(0)); + break; + } + case OpType::CX: { + apply_S(qbs.at(1)); + apply_V(qbs.at(1)); + apply_S(qbs.at(1)); + apply_CZ(qbs.at(0), qbs.at(1)); + apply_S(qbs.at(1)); + apply_V(qbs.at(1)); + apply_S(qbs.at(1)); + break; + } + case OpType::CY: { + apply_V(qbs.at(1)); + apply_CZ(qbs.at(0), qbs.at(1)); + apply_V(qbs.at(1)); + apply_V(qbs.at(1)); + apply_V(qbs.at(1)); + phase += 1.; + break; + } + case OpType::CZ: { + apply_CZ(qbs.at(0), qbs.at(1)); + break; + } + case OpType::ZZMax: { + apply_S(qbs.at(0)); + apply_S(qbs.at(1)); + apply_CZ(qbs.at(0), qbs.at(1)); + phase -= .25; + break; + } + case OpType::ECR: { + apply_S(qbs.at(0)); + apply_V(qbs.at(0)); + apply_V(qbs.at(0)); + apply_S(qbs.at(1)); + apply_V(qbs.at(1)); + apply_CZ(qbs.at(0), qbs.at(1)); + apply_S(qbs.at(1)); + apply_V(qbs.at(1)); + apply_S(qbs.at(1)); + phase += .25; + break; + } + case OpType::ISWAPMax: { + apply_V(qbs.at(0)); + apply_S(qbs.at(1)); + apply_V(qbs.at(1)); + apply_CZ(qbs.at(0), qbs.at(1)); + apply_V(qbs.at(0)); + apply_V(qbs.at(1)); + apply_CZ(qbs.at(0), qbs.at(1)); + apply_V(qbs.at(0)); + apply_V(qbs.at(1)); + apply_S(qbs.at(1)); + phase += 1.; + break; + } + case OpType::SWAP: { + apply_S(qbs.at(0)); + apply_V(qbs.at(0)); + apply_S(qbs.at(0)); + apply_CZ(qbs.at(0), qbs.at(1)); + apply_S(qbs.at(0)); + apply_V(qbs.at(0)); + apply_S(qbs.at(0)); + apply_S(qbs.at(1)); + apply_V(qbs.at(1)); + apply_S(qbs.at(1)); + apply_CZ(qbs.at(0), qbs.at(1)); + apply_S(qbs.at(0)); + apply_V(qbs.at(0)); + apply_S(qbs.at(0)); + apply_S(qbs.at(1)); + apply_V(qbs.at(1)); + apply_S(qbs.at(1)); + apply_CZ(qbs.at(0), qbs.at(1)); + apply_S(qbs.at(0)); + apply_V(qbs.at(0)); + apply_S(qbs.at(0)); + break; + } + case OpType::BRIDGE: { + apply_S(qbs.at(2)); + apply_V(qbs.at(2)); + apply_S(qbs.at(2)); + apply_CZ(qbs.at(0), qbs.at(2)); + apply_S(qbs.at(2)); + apply_V(qbs.at(2)); + apply_S(qbs.at(2)); + break; + } + case OpType::noop: { + break; + } + case OpType::Phase: { + throw std::logic_error( + "OpType::Phase cannot be applied via APState::apply_gate"); + } + default: { + throw BadOpType( + "Cannot be applied to a APState: not a Clifford gate", type); + } + } +} + +} // namespace tket \ No newline at end of file diff --git a/tket/src/Converters/APStateConverters.cpp b/tket/src/Converters/APStateConverters.cpp new file mode 100644 index 0000000000..a56e96da79 --- /dev/null +++ b/tket/src/Converters/APStateConverters.cpp @@ -0,0 +1,176 @@ +// Copyright 2019-2024 Cambridge Quantum Computing +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include "tket/Converters/Converters.hpp" + +namespace tket { + +APState circuit_to_apstate(const Circuit& circ) { + APState aps(circ.n_qubits()); + std::map qb_ordering; + for (const Qubit& q : circ.all_qubits()) + qb_ordering.insert({q, qb_ordering.size()}); + for (const Command& com : circ) { + auto args = com.get_args(); + std::vector qbs; + for (const UnitID& q : args) qbs.push_back(qb_ordering.at(q)); + aps.apply_gate(com.get_op_ptr()->get_type(), qbs); + } + return aps; +} + +Circuit apstate_to_circuit(const APState& ap) { + unsigned n_qbs = ap.A.rows(); + Circuit circ(n_qbs); + circ.qubit_create_all(); + for (unsigned q = 0; q < n_qbs; ++q) { + if (!ap.A(q, q)) + circ.add_op(OpType::H, {q}); + else if (ap.B(q)) + circ.add_op(OpType::X, {q}); + } + for (unsigned trgt = 0; trgt < n_qbs; ++trgt) { + if (ap.A(trgt, trgt)) { + for (unsigned ctrl = trgt + 1; ctrl < n_qbs; ++ctrl) + if (ap.A(trgt, ctrl)) circ.add_op(OpType::CX, {ctrl, trgt}); + } + } + for (unsigned q1 = 0; q1 < n_qbs; ++q1) { + if (ap.A(q1, q1)) continue; + for (unsigned q2 = q1 + 1; q2 < n_qbs; ++q2) { + if (ap.E(q1, q2)) circ.add_op(OpType::CZ, {q1, q2}); + } + switch (ap.P(q1) % 4) { + case 1: { + circ.add_op(OpType::S, {q1}); + break; + } + case 2: { + circ.add_op(OpType::Z, {q1}); + break; + } + case 3: { + circ.add_op(OpType::Sdg, {q1}); + break; + } + default: { + break; + } + } + } + circ.add_phase(ap.phase); + return circ; +} + +APState tableau_to_apstate(SymplecticTableau tab) { + unsigned n_qbs = tab.get_n_qubits(); + if (tab.get_n_rows() != n_qbs) + throw std::logic_error( + "tableau_to_apstate requires a tableau with n commuting rows for n " + "qubits"); + MatrixXb fullmat = MatrixXb::Zero(n_qbs, 2 * n_qbs); + /** + * Gaussian elimination by the x matrix first ensures the bottom rows are only + * Zs, i.e. describing rows of A. Reversing the columns of the x matrix + * guarantees that each row has an X on at most one free qubit, simplifying + * the code for finding E and P + */ + for (unsigned c = 0; c < n_qbs; ++c) { + fullmat.col(c) = tab.xmat.col(n_qbs - 1 - c); + } + fullmat.block(0, n_qbs, n_qbs, n_qbs) = tab.zmat; + std::vector> row_ops = + gaussian_elimination_row_ops(fullmat); + for (const std::pair& op : row_ops) + tab.row_mult(op.first, op.second); + + MatrixXb A = MatrixXb::Zero(n_qbs, n_qbs); + VectorXb B = VectorXb::Zero(n_qbs); + MatrixXb E = MatrixXb::Zero(n_qbs, n_qbs); + Eigen::VectorXi P = Eigen::VectorXi::Zero(n_qbs); + + unsigned n_leading = 0; + for (unsigned r = n_qbs; r > 0;) { + --r; + bool is_a_row = true; + for (unsigned c = 0; c < n_qbs; ++c) { + if (tab.xmat(r, c)) { + is_a_row = false; + break; + } + } + if (!is_a_row) break; + unsigned first = 0; + for (unsigned c = 0; c < n_qbs; ++c) { + if (tab.zmat(r, c)) { + first = c; + break; + } + } + A.row(first) = tab.zmat.row(r); + B(first) = tab.phase(r); + ++n_leading; + } + /** + * Each free qubit is after all leaders connected to it, by reduced + * row-echelon of A. Therefore Gaussian elimination of the x matrix in reverse + * order gives rows corresponding to the columns of A of free qubits (plus the + * free qubit itself). + * + * Then the corresponding row of the z matrix is exactly the free qubit's row + * of E, from which we take the diagonal and phase vector to determine P. + */ + for (unsigned r = 0; r < n_qbs - n_leading; ++r) { + unsigned free = 0; + for (unsigned c = n_qbs; c > 0;) { + --c; + if (tab.xmat(r, c)) { + free = c; + break; + } + } + E.row(free) = tab.zmat.row(r); + if (tab.zmat(r, free)) P(free) += 1; + if (tab.phase(r)) P(free) += 2; + E(free, free) = false; + } + + return APState(A, B, E, P, 0); +} + +SymplecticTableau apstate_to_tableau(const APState& ap) { + unsigned n_qbs = ap.A.rows(); + MatrixXb xmat = MatrixXb::Zero(n_qbs, n_qbs); + MatrixXb zmat = MatrixXb::Zero(n_qbs, n_qbs); + VectorXb phase = VectorXb::Zero(n_qbs); + + for (unsigned q = 0; q < n_qbs; ++q) { + if (ap.A(q, q)) { + // Leading qubit + zmat.row(q) = ap.A.row(q); + phase(q) = ap.B(q); + } else { + // Free qubit + xmat.row(q) = ap.A.col(q); + xmat(q, q) = true; + zmat.row(q) = ap.E.row(q); + zmat(q, q) = (ap.P(q) % 2 == 1); + phase(q) = (ap.P(q) % 4 > 1); + } + } + + return SymplecticTableau(xmat, zmat, phase); +} + +} // namespace tket \ No newline at end of file diff --git a/tket/test/CMakeLists.txt b/tket/test/CMakeLists.txt index 4e020c4d16..fe7011d6c5 100644 --- a/tket/test/CMakeLists.txt +++ b/tket/test/CMakeLists.txt @@ -102,6 +102,7 @@ add_executable(test-tket src/Simulation/test_CircuitSimulator.cpp src/Simulation/test_PauliExpBoxUnitaryCalculator.cpp src/test_AASRoute.cpp + src/test_APState.cpp src/test_ArchitectureAwareSynthesis.cpp src/test_Architectures.cpp src/test_Assertion.cpp diff --git a/tket/test/src/test_APState.cpp b/tket/test/src/test_APState.cpp new file mode 100644 index 0000000000..12320e8a90 --- /dev/null +++ b/tket/test/src/test_APState.cpp @@ -0,0 +1,347 @@ +// Copyright 2019-2024 Cambridge Quantum Computing +// +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include + +#include "Simulation/ComparisonFunctions.hpp" +#include "testutil.hpp" +#include "tket/Circuit/Simulation/CircuitSimulator.hpp" +#include "tket/Clifford/APState.hpp" +#include "tket/Converters/Converters.hpp" + +namespace tket { +namespace test_APState { + +void test_apply_gate( + const MatrixXb& A, const VectorXb& B, const MatrixXb& E, + const Eigen::VectorXi& P, OpType ot, const std::vector& args) { + APState ap(A, B, E, P, 0); + ap.verify(); + auto sv_before = ap.to_statevector(); + + Circuit circ(A.rows()); + circ.add_op(ot, args); + auto gate_u = tket_sim::get_unitary(circ); + + ap.apply_gate(ot, args); + + ap.verify(); + auto sv_after = ap.to_statevector(); + + CHECK(tket_sim::compare_statevectors_or_unitaries( + gate_u * sv_before, sv_after, tket_sim::MatrixEquivalence::EQUAL)); +} + +SCENARIO("CZ cases") { + GIVEN("CZ on free qubits") { + MatrixXb A = MatrixXb::Zero(4, 4); + VectorXb B = VectorXb::Zero(4); + MatrixXb E = MatrixXb::Zero(4, 4); + Eigen::VectorXi P = Eigen::VectorXi::Zero(4); + A(0, 0) = A(0, 1) = A(0, 3) = true; + E(2, 3) = E(3, 2) = true; + test_apply_gate(A, B, E, P, OpType::CZ, {1, 2}); + } + GIVEN("CZ on one leading qubit and connected free") { + for (unsigned b = 0; b < 2; ++b) { + MatrixXb A = MatrixXb::Zero(5, 5); + VectorXb B = VectorXb::Zero(5); + MatrixXb E = MatrixXb::Zero(5, 5); + Eigen::VectorXi P = Eigen::VectorXi::Zero(5); + A(0, 0) = A(0, 2) = A(0, 3) = true; + A(1, 1) = A(1, 3) = A(1, 4) = true; + B(1) = (b == 1); + test_apply_gate(A, B, E, P, OpType::CZ, {0, 3}); + } + } + GIVEN("CZ on one leading qubit and unconnected free") { + for (unsigned b = 0; b < 2; ++b) { + MatrixXb A = MatrixXb::Zero(5, 5); + VectorXb B = VectorXb::Zero(5); + MatrixXb E = MatrixXb::Zero(5, 5); + Eigen::VectorXi P = Eigen::VectorXi::Zero(5); + A(0, 0) = A(0, 2) = A(0, 3) = true; + A(1, 1) = A(1, 3) = A(1, 4) = true; + B(1) = (b == 1); + test_apply_gate(A, B, E, P, OpType::CZ, {0, 4}); + } + } + GIVEN("CZ on leading qubits") { + for (unsigned b1 = 0; b1 < 2; ++b1) { + for (unsigned b2 = 0; b2 < 2; ++b2) { + MatrixXb A = MatrixXb::Zero(8, 8); + VectorXb B = VectorXb::Zero(8); + MatrixXb E = MatrixXb::Zero(8, 8); + Eigen::VectorXi P = Eigen::VectorXi::Zero(8); + A(0, 0) = A(0, 2) = A(0, 3) = A(0, 4) = A(0, 5) = true; + A(1, 1) = A(1, 4) = A(1, 5) = A(1, 6) = A(1, 7) = true; + B(0) = (b1 == 1); + B(1) = (b2 == 1); + test_apply_gate(A, B, E, P, OpType::CZ, {0, 1}); + } + } + } +} + +SCENARIO("S cases") { + GIVEN("S on free qubit") { + MatrixXb A = MatrixXb::Zero(3, 3); + VectorXb B = VectorXb::Zero(3); + MatrixXb E = MatrixXb::Zero(3, 3); + Eigen::VectorXi P = Eigen::VectorXi::Zero(3); + A(0, 0) = A(0, 1) = A(0, 2) = true; + E(1, 2) = E(2, 1) = true; + test_apply_gate(A, B, E, P, OpType::S, {2}); + } + GIVEN("S on leading qubit") { + for (unsigned b = 0; b < 2; ++b) { + MatrixXb A = MatrixXb::Zero(4, 4); + VectorXb B = VectorXb::Zero(4); + MatrixXb E = MatrixXb::Zero(4, 4); + Eigen::VectorXi P = Eigen::VectorXi::Zero(4); + A(0, 0) = A(0, 1) = A(0, 2) = true; + B(0) = (b == 1); + E(1, 3) = E(3, 1) = true; + test_apply_gate(A, B, E, P, OpType::S, {0}); + } + } + GIVEN("S on disconnected leading qubit") { + for (unsigned b = 0; b < 2; ++b) { + MatrixXb A = MatrixXb::Zero(1, 1); + VectorXb B = VectorXb::Zero(1); + MatrixXb E = MatrixXb::Zero(1, 1); + Eigen::VectorXi P = Eigen::VectorXi::Zero(1); + A(0, 0) = true; + B(0) = (b == 1); + test_apply_gate(A, B, E, P, OpType::S, {0}); + } + } +} + +SCENARIO("V cases") { + GIVEN("V on leading qubit") { + for (unsigned b = 0; b < 2; ++b) { + MatrixXb A = MatrixXb::Zero(4, 4); + VectorXb B = VectorXb::Zero(4); + MatrixXb E = MatrixXb::Zero(4, 4); + Eigen::VectorXi P = Eigen::VectorXi::Zero(4); + A(0, 0) = A(0, 2) = A(0, 3) = true; + A(1, 1) = A(1, 3) = true; + B(0) = (b == 1); + test_apply_gate(A, B, E, P, OpType::V, {0}); + } + } + GIVEN("V on free qubit with some leading") { + for (unsigned b = 0; b < 2; ++b) { + for (unsigned p = 0; p < 4; ++p) { + MatrixXb A = MatrixXb::Zero(9, 9); + VectorXb B = VectorXb::Zero(9); + MatrixXb E = MatrixXb::Zero(9, 9); + Eigen::VectorXi P = Eigen::VectorXi::Zero(9); + A(0, 0) = A(0, 2) = A(0, 4) = A(0, 5) = A(0, 7) = true; + A(1, 1) = A(1, 2) = A(1, 3) = A(1, 4) = A(1, 5) = A(1, 6) = true; + B(1) = (b == 1); + E(4, 5) = E(5, 4) = true; + E(4, 6) = E(6, 4) = true; + E(4, 7) = E(7, 4) = true; + E(4, 8) = E(8, 4) = true; + P(4) = p; + test_apply_gate(A, B, E, P, OpType::V, {4}); + } + } + } + GIVEN("V on free qubit with some earlier connected free") { + for (unsigned p1 = 0; p1 < 4; ++p1) { + for (unsigned p2 = 0; p2 < 4; ++p2) { + MatrixXb A = MatrixXb::Zero(9, 9); + VectorXb B = VectorXb::Zero(9); + MatrixXb E = MatrixXb::Zero(9, 9); + Eigen::VectorXi P = Eigen::VectorXi::Zero(9); + A(0, 0) = true; + A(1, 1) = A(1, 4) = A(1, 6) = true; + A(2, 2) = A(2, 4) = A(2, 7) = true; + A(3, 3) = A(3, 4) = A(3, 8) = A(3, 8) = true; + E(4, 5) = E(5, 4) = true; + E(4, 7) = E(7, 4) = true; + E(4, 8) = E(8, 4) = true; + E(5, 6) = E(6, 5) = true; + E(5, 7) = E(7, 5) = true; + E(5, 8) = E(8, 5) = true; + P(4) = p1; + P(5) = p2; + test_apply_gate(A, B, E, P, OpType::V, {5}); + } + } + } + GIVEN("V on free qubit with no earlier connected free") { + for (unsigned p = 0; p < 4; ++p) { + MatrixXb A = MatrixXb::Zero(4, 4); + VectorXb B = VectorXb::Zero(4); + MatrixXb E = MatrixXb::Zero(4, 4); + Eigen::VectorXi P = Eigen::VectorXi::Zero(4); + A(0, 0) = A(0, 2) = true; + E(1, 2) = E(2, 1) = true; + E(1, 3) = E(3, 1) = true; + P(1) = p; + test_apply_gate(A, B, E, P, OpType::V, {1}); + } + } + GIVEN("V on disconnected free qubit") { + for (unsigned p = 0; p < 4; ++p) { + MatrixXb A = MatrixXb::Zero(1, 1); + VectorXb B = VectorXb::Zero(1); + MatrixXb E = MatrixXb::Zero(1, 1); + Eigen::VectorXi P = Eigen::VectorXi::Zero(1, 1); + P(0) = p; + test_apply_gate(A, B, E, P, OpType::V, {0}); + } + } +} + +SCENARIO("Gate encodings") { + std::list>> test_gates = { + {OpType::Z, {0}}, {OpType::X, {0}}, + {OpType::Y, {0}}, {OpType::S, {0}}, + {OpType::Sdg, {0}}, {OpType::V, {0}}, + {OpType::Vdg, {0}}, {OpType::SX, {0}}, + {OpType::SXdg, {0}}, {OpType::H, {0}}, + {OpType::CX, {0, 1}}, {OpType::CY, {0, 1}}, + {OpType::CZ, {0, 1}}, {OpType::ZZMax, {0, 1}}, + {OpType::ECR, {0, 1}}, {OpType::ISWAPMax, {0, 1}}, + {OpType::SWAP, {0, 1}}, {OpType::BRIDGE, {0, 1, 2}}, + {OpType::noop, {0}}, + }; + GIVEN("Check Z actions") { + for (const auto& com : test_gates) { + MatrixXb A = MatrixXb::Identity(3, 3); + VectorXb B = VectorXb::Zero(3); + MatrixXb E = MatrixXb::Zero(3, 3); + Eigen::VectorXi P = Eigen::VectorXi::Zero(3); + test_apply_gate(A, B, E, P, com.first, com.second); + } + } + GIVEN("Check X actions") { + for (const auto& com : test_gates) { + MatrixXb A = MatrixXb::Zero(3, 3); + VectorXb B = VectorXb::Zero(3); + MatrixXb E = MatrixXb::Zero(3, 3); + Eigen::VectorXi P = Eigen::VectorXi::Zero(3); + test_apply_gate(A, B, E, P, com.first, com.second); + } + } +} + +SCENARIO("Loading from a statevector") { + MatrixXb A = MatrixXb::Zero(4, 4); + VectorXb B = VectorXb::Zero(4); + MatrixXb E = MatrixXb::Zero(4, 4); + Eigen::VectorXi P = Eigen::VectorXi::Zero(4); + A(0, 0) = A(0, 2) = A(0, 3) = true; + A(1, 1) = A(1, 2) = true; + B(0) = true; + E(2, 3) = E(3, 2) = true; + P(2) = 1; + P(3) = 2; + APState ap(A, B, E, P, 0); + auto sv = ap.to_statevector(); + APState reconstructed(sv); + auto sv2 = reconstructed.to_statevector(); + CHECK(tket_sim::compare_statevectors_or_unitaries( + sv, sv2, tket_sim::MatrixEquivalence::EQUAL)); +} + +SCENARIO("Converting from/to a circuit") { + GIVEN("A circuit in the standard AP form") { + Circuit circ(4); + circ.qubit_create_all(); + circ.add_op(OpType::H, {2}); + circ.add_op(OpType::H, {3}); + circ.add_op(OpType::CX, {2, 0}); + circ.add_op(OpType::CX, {2, 1}); + circ.add_op(OpType::CX, {3, 1}); + circ.add_op(OpType::S, {2}); + circ.add_op(OpType::Z, {3}); + APState ap = circuit_to_apstate(circ); + auto sv_circ = tket_sim::get_statevector(circ); + auto sv_ap = ap.to_statevector(); + CHECK(tket_sim::compare_statevectors_or_unitaries( + sv_circ, sv_ap, tket_sim::MatrixEquivalence::EQUAL)); + Circuit reconstructed = apstate_to_circuit(ap); + CHECK(circ == reconstructed); + } + GIVEN("A generic circuit") { + Circuit circ(4); + circ.qubit_create_all(); + circ.add_op(OpType::V, {0}); + circ.add_op(OpType::CX, {0, 1}); + circ.add_op(OpType::CY, {1, 3}); + circ.add_op(OpType::H, {3}); + circ.add_op(OpType::ZZMax, {2, 3}); + APState ap = circuit_to_apstate(circ); + auto sv_circ = tket_sim::get_statevector(circ); + auto sv_ap = ap.to_statevector(); + CHECK(tket_sim::compare_statevectors_or_unitaries( + sv_circ, sv_ap, tket_sim::MatrixEquivalence::EQUAL)); + Circuit reconstructed = apstate_to_circuit(ap); + auto sv_rec = tket_sim::get_statevector(reconstructed); + CHECK(tket_sim::compare_statevectors_or_unitaries( + sv_circ, sv_rec, tket_sim::MatrixEquivalence::EQUAL)); + } +} + +SCENARIO("Converting from/to a tableau") { + GIVEN("Check up to global phase using circuit") { + Circuit circ(8); + circ.qubit_create_all(); + circ.add_op(OpType::X, {1}); + circ.add_op(OpType::X, {5}); + circ.add_op(OpType::H, {2}); + circ.add_op(OpType::H, {4}); + circ.add_op(OpType::H, {6}); + circ.add_op(OpType::H, {7}); + circ.add_op(OpType::CX, {2, 1}); + circ.add_op(OpType::CX, {4, 0}); + circ.add_op(OpType::CX, {4, 3}); + circ.add_op(OpType::CX, {6, 0}); + circ.add_op(OpType::CX, {6, 1}); + circ.add_op(OpType::CX, {7, 5}); + circ.add_op(OpType::CZ, {2, 6}); + circ.add_op(OpType::CZ, {4, 6}); + circ.add_op(OpType::CZ, {4, 7}); + circ.add_op(OpType::CZ, {6, 7}); + circ.add_op(OpType::S, {2}); + circ.add_op(OpType::Sdg, {4}); + circ.add_op(OpType::Z, {7}); + ChoiMixTableau cmt = circuit_to_cm_tableau(circ); + APState ap = tableau_to_apstate(cmt.tab_); + auto sv_circ = tket_sim::get_statevector(circ); + auto sv_ap = ap.to_statevector(); + CHECK(tket_sim::compare_statevectors_or_unitaries( + sv_circ, sv_ap, tket_sim::MatrixEquivalence::EQUAL)); + SymplecticTableau tab2 = apstate_to_tableau(ap); + ChoiMixTableau cmt2(tab2.xmat, tab2.zmat, tab2.phase); + auto [circ2, perm] = cm_tableau_to_exact_circuit(cmt2); + qubit_map_t inv; + for (const std::pair& qp : perm) + inv.insert({qp.second, qp.first}); + circ2.permute_boundary_output(inv); + auto sv_circ2 = tket_sim::get_statevector(circ2); + CHECK(tket_sim::compare_statevectors_or_unitaries( + sv_circ, sv_circ2, + tket_sim::MatrixEquivalence::EQUAL_UP_TO_GLOBAL_PHASE)); + } +} + +} // namespace test_APState +} // namespace tket \ No newline at end of file From 8dfb97a319c4eb37f17dabbf471f8f57c7c6e86e Mon Sep 17 00:00:00 2001 From: Will Simmons Date: Wed, 13 Nov 2024 09:12:34 +0000 Subject: [PATCH 2/6] Attempt signed shift fix --- tket/src/Clifford/APState.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tket/src/Clifford/APState.cpp b/tket/src/Clifford/APState.cpp index a21e4b446b..2dfa489611 100644 --- a/tket/src/Clifford/APState.cpp +++ b/tket/src/Clifford/APState.cpp @@ -75,7 +75,7 @@ APState::APState(const Eigen::VectorXcd& sv) { for (unsigned x = 1; x < sv.size(); ++x) { if (sv(z0 ^ x) != 0.) { ++n_non_zero; - if (n_non_zero == (1 << offsets.size())) offsets.push_back(x); + if (n_non_zero == (1u << offsets.size())) offsets.push_back(x); } } @@ -216,7 +216,7 @@ Eigen::VectorXcd APState::to_statevector() { Eigen::VectorXcd sv = Eigen::VectorXcd::Zero(1 << n_qubits); unsigned n_terms = 0; Complex g_phase = std::exp(i_ * PI * eval_expr(phase).value()); - for (unsigned x = 0; x < 1 << n_qubits; ++x) { + for (unsigned x = 0; x < (1u << n_qubits); ++x) { VectorXb x_binary = VectorXb::Zero(n_qubits); for (unsigned i = 0; i < n_qubits; ++i) { unsigned mask = 1 << (n_qubits - 1 - i); From e51bd2e09f2bdf90060a0b93efad18cd9a168ae1 Mon Sep 17 00:00:00 2001 From: Will Simmons Date: Wed, 13 Nov 2024 09:12:44 +0000 Subject: [PATCH 3/6] Bump tket version --- tket/conanfile.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tket/conanfile.py b/tket/conanfile.py index 147e674cb2..31633eef95 100644 --- a/tket/conanfile.py +++ b/tket/conanfile.py @@ -23,7 +23,7 @@ class TketConan(ConanFile): name = "tket" - version = "1.3.39" + version = "1.3.40" package_type = "library" license = "Apache 2" homepage = "https://github.com/CQCL/tket" From dc60ac30a8230546bfd7661c445b8fde9b51dad8 Mon Sep 17 00:00:00 2001 From: Will Simmons Date: Wed, 13 Nov 2024 09:21:56 +0000 Subject: [PATCH 4/6] Forgot to change the version number in pytket, oops! --- pytket/conanfile.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pytket/conanfile.py b/pytket/conanfile.py index e8b79f648b..9b76b1cb1a 100644 --- a/pytket/conanfile.py +++ b/pytket/conanfile.py @@ -38,7 +38,7 @@ def requirements(self): self.requires("pybind11_json/0.2.14") self.requires("symengine/0.13.0") self.requires("tkassert/0.3.4@tket/stable") - self.requires("tket/1.3.46@tket/stable") + self.requires("tket/1.3.47@tket/stable") self.requires("tklog/0.3.3@tket/stable") self.requires("tkrng/0.3.3@tket/stable") self.requires("tktokenswap/0.3.9@tket/stable") From f1881a743614c3a5495ff57a02f480d5ba46f5d7 Mon Sep 17 00:00:00 2001 From: Will Simmons Date: Thu, 14 Nov 2024 09:12:22 +0000 Subject: [PATCH 5/6] Doxygen and type warnings fix --- tket/include/tket/Clifford/APState.hpp | 2 +- tket/src/Converters/APStateConverters.cpp | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/tket/include/tket/Clifford/APState.hpp b/tket/include/tket/Clifford/APState.hpp index 402106b11c..4d9577c4e8 100644 --- a/tket/include/tket/Clifford/APState.hpp +++ b/tket/include/tket/Clifford/APState.hpp @@ -36,7 +36,7 @@ namespace tket { * CZ gates and a (n-k) vector P of integers mod 4 describing S gates. * * This gives a canonical statevector: - * \sum_{x, Ax=B} i^{\Phi(x)} \ket{x} + * \sum_{x, Ax=B} i^{\Phi(x)} |{x}> * * This canonical statevector fixes a reference phase from which we can track * global phase with an additional parameter. diff --git a/tket/src/Converters/APStateConverters.cpp b/tket/src/Converters/APStateConverters.cpp index a56e96da79..d73dcbc53d 100644 --- a/tket/src/Converters/APStateConverters.cpp +++ b/tket/src/Converters/APStateConverters.cpp @@ -20,7 +20,7 @@ APState circuit_to_apstate(const Circuit& circ) { APState aps(circ.n_qubits()); std::map qb_ordering; for (const Qubit& q : circ.all_qubits()) - qb_ordering.insert({q, qb_ordering.size()}); + qb_ordering.insert({q, (unsigned)qb_ordering.size()}); for (const Command& com : circ) { auto args = com.get_args(); std::vector qbs; From 8d7dc861c3c61bb1bfdf6ea4da41a413914d5623 Mon Sep 17 00:00:00 2001 From: Will Simmons Date: Thu, 14 Nov 2024 09:40:24 +0000 Subject: [PATCH 6/6] Missed a second doxygen error --- tket/include/tket/Clifford/APState.hpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tket/include/tket/Clifford/APState.hpp b/tket/include/tket/Clifford/APState.hpp index 4d9577c4e8..63b90ecef2 100644 --- a/tket/include/tket/Clifford/APState.hpp +++ b/tket/include/tket/Clifford/APState.hpp @@ -36,7 +36,7 @@ namespace tket { * CZ gates and a (n-k) vector P of integers mod 4 describing S gates. * * This gives a canonical statevector: - * \sum_{x, Ax=B} i^{\Phi(x)} |{x}> + * \sum_{x, Ax=B} i^{\Phi(x)} |x> * * This canonical statevector fixes a reference phase from which we can track * global phase with an additional parameter. @@ -101,7 +101,7 @@ class APState { const Eigen::VectorXi& P_, const Expr& phase_); /** - * Constructs the state \ket{0}^{\otimes n_qubits} in AP form. + * Constructs the state |0>^{\otimes n_qubits} in AP form. */ APState(unsigned n_qubits);