diff --git a/CHANGELOG.md b/CHANGELOG.md index ca9380e82..747b234e4 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -11,6 +11,14 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - #476: Introduced new methods `OpenGraph.extract_circuit`, `CliffordMap.to_tableau` and new function `graphix.circ_ext.compilation.cm_berg_pass`. Circuit extraction can be done natively in Graphix. +- #484: J & CZ transpilation. + - Replaced `Circuit.transpile()` with a new approach based decomposing circuits into J & CZ gates. + - Added `Circuit.transpile_to_causal_flow()` to produce `CausalFlow` using the same decomposition. + - Added `instruction.J` class. + - Added `Circuit.transpile_j_to_rzh()` method to prepare circuits with J gates for export to OpenQASM. + - Added `transpile` argument to `qasm3_exporter.circuit_to_qasm3` and `qasm3_exporter.circuit_to_qasm3_lines`, which defaults to true and applies `Circuit.transpile_j_to_rzh` and `Circuit.transpile_measurements_to_z_axis` methods. + - The transpiler now returns `TranspiledPattern` or `TranspiledFlow`, instead of `TranspileResult`. + - #490: Introduced new `Instruction` and `Command` namespace classes for instruction and command instantiation. - #505 diff --git a/docs/source/generator.rst b/docs/source/generator.rst index 1a9915952..97918ae39 100644 --- a/docs/source/generator.rst +++ b/docs/source/generator.rst @@ -16,6 +16,8 @@ Pattern Generation .. automethod:: simulate_statevector + .. automethod:: cz + .. automethod:: cnot .. automethod:: h @@ -34,10 +36,14 @@ Pattern Generation .. automethod:: rz + .. automethod:: j + .. automethod:: ccx .. automethod:: m -.. autoclass:: TranspileResult +.. autoclass:: TranspiledPattern + +.. autoclass:: TranspiledFlow .. autoclass:: SimulateResult diff --git a/examples/deutsch_jozsa.py b/examples/deutsch_jozsa.py index 68a032846..5f4d6dee4 100644 --- a/examples/deutsch_jozsa.py +++ b/examples/deutsch_jozsa.py @@ -73,6 +73,7 @@ # %% # Now we preprocess all Pauli measurements, which requires that we move inputs to N commands +pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() print( pattern.to_ascii( diff --git a/examples/ghz_with_tn.py b/examples/ghz_with_tn.py index 3e73e8364..2834c4e0e 100644 --- a/examples/ghz_with_tn.py +++ b/examples/ghz_with_tn.py @@ -16,7 +16,7 @@ from graphix import Circuit -n = 100 +n = 50 print(f"{n}-qubit GHZ state generation") circuit = Circuit(n) diff --git a/examples/mbqc_vqe.py b/examples/mbqc_vqe.py index 9b690df9f..6f991bb20 100644 --- a/examples/mbqc_vqe.py +++ b/examples/mbqc_vqe.py @@ -92,7 +92,9 @@ def build_mbqc_pattern(self, params: Iterable[ParameterizedAngle]) -> Pattern: pattern = circuit.transpile().pattern pattern.standardize() pattern.shift_signals() - pattern.remove_pauli_measurements() # Perform Pauli measurements + pattern = pattern.infer_pauli_measurements() + pattern.remove_pauli_measurements() + pattern.minimize_space() return pattern # %% @@ -155,7 +157,7 @@ def cost_function(params: Iterable[float]) -> float: # %% # Perform the optimization using COBYLA def compute() -> OptimizeResult: - return minimize(cost_function, initial_params, method="COBYLA", options={"maxiter": 100}) + return minimize(cost_function, initial_params, method="COBYLA", options={"maxiter": 20}) result = compute() @@ -172,9 +174,9 @@ def compute() -> OptimizeResult: # Compare performances between using parameterized circuits (with placeholders) or not mbqc_vqe = MBQCVQEWithPlaceholders(n_qubits, hamiltonian) -time_with_placeholders = timeit(compute, number=2) +time_with_placeholders = timeit(compute, number=1) print(f"Time with placeholders: {time_with_placeholders}") mbqc_vqe = MBQCVQE(n_qubits, hamiltonian) -time_without_placeholders = timeit(compute, number=2) +time_without_placeholders = timeit(compute, number=1) print(f"Time without placeholders: {time_without_placeholders}") diff --git a/examples/qaoa.py b/examples/qaoa.py index f0f2c1a42..da993fe37 100644 --- a/examples/qaoa.py +++ b/examples/qaoa.py @@ -40,6 +40,7 @@ # %% # perform Pauli measurements and plot the new (minimal) graph to perform the same quantum computation +pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() pattern.draw(flow_from_pattern=False) diff --git a/examples/qft_with_tn.py b/examples/qft_with_tn.py index 8abb6ded5..5cb991c2a 100644 --- a/examples/qft_with_tn.py +++ b/examples/qft_with_tn.py @@ -67,6 +67,7 @@ def qft(circuit: Circuit, n: int) -> None: # Using graph rewriting rules, we can classically preprocess Pauli measurements. # We are currently improving the speed of this process by using rust-based graph manipulation backend. pattern.remove_input_nodes() +pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements(standardize=True) diff --git a/examples/tn_simulation.py b/examples/tn_simulation.py index 63c061e1b..84a5e8b19 100644 --- a/examples/tn_simulation.py +++ b/examples/tn_simulation.py @@ -85,6 +85,7 @@ def ansatz( # %% # Optimizing by removing Pauli measurements in the pattern. pattern.remove_input_nodes() +pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements(standardize=True) # %% @@ -203,6 +204,7 @@ def cost( pattern.standardize() pattern.shift_signals() pattern.remove_input_nodes() + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements(standardize=True) mbqc_tn = pattern.simulate_pattern(backend="tensornetwork", graph_prep="parallel") exp_val: float = 0 diff --git a/examples/visualization.py b/examples/visualization.py index 8b4e8742e..d6e0e68de 100644 --- a/examples/visualization.py +++ b/examples/visualization.py @@ -37,6 +37,7 @@ # %% # next, show the gflow: +pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() pattern.draw(flow_from_pattern=False, measurement_labels=True) diff --git a/graphix/instruction.py b/graphix/instruction.py index f608ccf8f..86be7e347 100644 --- a/graphix/instruction.py +++ b/graphix/instruction.py @@ -49,6 +49,7 @@ class InstructionKind(Enum): X = enum.auto() Y = enum.auto() Z = enum.auto() + J = enum.auto() I = enum.auto() M = enum.auto() RX = enum.auto() @@ -290,6 +291,19 @@ def visit(self, visitor: InstructionVisitor) -> RZ: return RZ(visitor.visit_qubit(self.target), visitor.visit_angle(self.angle)) +@dataclass(repr=False) +class J(_KindChecker, BaseInstruction): + """J circuit instruction.""" + + target: int + angle: ParameterizedAngle = field(metadata={"repr": repr_angle}) + kind: ClassVar[Literal[InstructionKind.J]] = field(default=InstructionKind.J, init=False) + + @override + def visit(self, visitor: InstructionVisitor) -> J: + return J(visitor.visit_qubit(self.target), visitor.visit_angle(self.angle)) + + class InstructionWithoutRZZ: """Grouping of all instructions except RZZ for namespace exposure. @@ -313,6 +327,7 @@ class InstructionWithoutRZZ: RX: TypeAlias = RX RY: TypeAlias = RY RZ: TypeAlias = RZ + J: TypeAlias = J def __init__(self) -> None: raise TypeError("InstructionWithoutRZZ is a namespace, not a class.") @@ -334,5 +349,5 @@ def __init__(self) -> None: if TYPE_CHECKING: - InstructionTypeWithoutRZZ = CCX | CNOT | SWAP | CZ | H | S | X | Y | Z | I | M | RX | RY | RZ + InstructionTypeWithoutRZZ = CCX | CNOT | SWAP | CZ | H | S | X | Y | Z | I | M | RX | RY | RZ | J InstructionType = InstructionTypeWithoutRZZ | RZZ diff --git a/graphix/ops.py b/graphix/ops.py index 075273da7..381a212df 100644 --- a/graphix/ops.py +++ b/graphix/ops.py @@ -167,6 +167,35 @@ def rz(theta: ParameterizedAngle) -> npt.NDArray[np.complex128] | npt.NDArray[np """ return Ops._cast_array([[exp(-1j * angle_to_rad(theta) / 2), 0], [0, exp(1j * angle_to_rad(theta) / 2)]], theta) + @overload + @staticmethod + def j(theta: Angle) -> npt.NDArray[np.complex128]: ... + + @overload + @staticmethod + def j(theta: Expression) -> npt.NDArray[np.object_]: ... + + @staticmethod + def j(theta: ParameterizedAngle) -> npt.NDArray[np.complex128] | npt.NDArray[np.object_]: + """J operation. + + Parameters + ---------- + theta : Angle | Expression + rotation angle in units of π + + Returns + ------- + operator : 2*2 np.asarray + """ + return Ops._cast_array( + [ + [1 / np.sqrt(2), (1 / np.sqrt(2)) * exp(1j * angle_to_rad(theta))], + [1 / np.sqrt(2), (-1 / np.sqrt(2)) * exp(1j * angle_to_rad(theta))], + ], + theta, + ) + @overload @staticmethod def rzz(theta: Angle) -> npt.NDArray[np.complex128]: ... diff --git a/graphix/pattern.py b/graphix/pattern.py index bcc18035b..e632a164b 100644 --- a/graphix/pattern.py +++ b/graphix/pattern.py @@ -388,16 +388,16 @@ def to_ascii( >>> circuit.ccx(0, 1, 2) >>> pattern = circuit.transpile().pattern >>> pattern.to_ascii() - 'Z(16,{11}) Z(18,{13}) Z(20,{7}) X(16,{15}) X(18,{14}) X(20,{19}) {17}[M(19)]{7} {10}[M(15,-pi/4)]{11} {17,12}[M(14)]{13} {17,9}[M(11)]{10} {0,1,5,10,13}[M(7,-pi/4)]{17} {1}[M(13,-7pi/4)]{12} {8}[M(10,-7pi/4)]{9} {17}[M(12)]{1} {6}[M(9)]{8} {8,3}[M(1,-pi/4)] {5}[M(8,-pi/4)]{6} {17,4}[M(6)]{5} [M(17)]{0} {3}[M(5,-7pi/4)]{4} M(0) {2}[M(4)]{3} [M(3)]{2} M(2) E(19,20) E(15,16) E(14,18) E(11,15) E(7,19) E(7,14) E(7,11) E(13,14) E(10,11) E(12,13) E(12,7) E(9,10) E(1,12) E(1,9) E(8,9)...(27 more commands)' + 'X(24,{21}) Z(24,{20}) X(32,{31}) Z(32,{28}) X(30,{29}) Z(30,{0}) {27}[M(31,0)]{28} E(31,32) N(32) [M(29,0)]{0} E(29,30) N(30) {26}[M(28,3pi/2)]{27} E(28,31) N(31) {17}[M(21,0)]{20} E(21,24) N(24) {25}[M(27,0)]{26} E(27,28) N(28) {26,19,15,7}[M(0,7pi/4)] E(0,27) E(0,29) N(29) {16}[M(20,0)]{17} E(20,21) N(21) {23}[M(26,0)]{25} E(26,27) N(27) {22}[M(25,0)]{23} E(25,26) N(26) {15}[M(17,7pi/4)]{16} E(17,20) N(20) {19}[M(23,pi/4)]{22} E(23,25)...(63 more commands)' >>> pattern.to_ascii(left_to_right=True) - 'N(3) N(4) N(5) N(6) N(7) N(8) N(9) N(10) N(11) N(12) N(13) N(14) N(15) N(16) N(17) N(18) N(19) N(20) E(2,3) E(3,4) E(4,5) E(4,1) E(0,17) E(5,6) E(17,7) E(6,8) E(6,7) E(8,9) E(1,9) E(1,12) E(9,10) E(12,7) E(12,13) E(10,11) E(13,14) E(7,11) E(7,14) E(7,19) E(11,15)...(27 more commands)' + 'N(3) E(2,3) M(2,0) N(4) E(3,4) [M(3,0)]{2} E(4,1) N(5) E(4,5) {2}[M(4,0)]{3} N(6) E(5,6) {3}[M(5,pi/4)]{4} N(7) E(6,7) {4}[M(6,0)]{5} N(8) E(7,8) {5}[M(7,0)]{6} E(8,0) N(9) E(8,9) {6}[M(8,0)]{7} N(10) E(9,10) {7}[M(9,7pi/4)]{8} N(11) E(10,11) {8}[M(10,0)]{9} N(12) E(11,12) {9}[M(11,0)]{10} N(18) E(1,18) E(1,12) {11,3}[M(1,pi/4)] N(13) E(12,13) {10}[M(12,0)]{11}...(63 more commands)' >>> pattern.to_ascii(limit=None) - 'Z(16,{11}) Z(18,{13}) Z(20,{7}) X(16,{15}) X(18,{14}) X(20,{19}) {17}[M(19)]{7} {10}[M(15,-pi/4)]{11} {17,12}[M(14)]{13} {17,9}[M(11)]{10} {0,1,5,10,13}[M(7,-pi/4)]{17} {1}[M(13,-7pi/4)]{12} {8}[M(10,-7pi/4)]{9} {17}[M(12)]{1} {6}[M(9)]{8} {8,3}[M(1,-pi/4)] {5}[M(8,-pi/4)]{6} {17,4}[M(6)]{5} [M(17)]{0} {3}[M(5,-7pi/4)]{4} M(0) {2}[M(4)]{3} [M(3)]{2} M(2) E(19,20) E(15,16) E(14,18) E(11,15) E(7,19) E(7,14) E(7,11) E(13,14) E(10,11) E(12,13) E(12,7) E(9,10) E(1,12) E(1,9) E(8,9) E(6,7) E(6,8) E(17,7) E(5,6) E(0,17) E(4,1) E(4,5) E(3,4) E(2,3) N(20) N(19) N(18) N(17) N(16) N(15) N(14) N(13) N(12) N(11) N(10) N(9) N(8) N(7) N(6) N(5) N(4) N(3)' + 'X(24,{21}) Z(24,{20}) X(32,{31}) Z(32,{28}) X(30,{29}) Z(30,{0}) {27}[M(31,0)]{28} E(31,32) N(32) [M(29,0)]{0} E(29,30) N(30) {26}[M(28,3pi/2)]{27} E(28,31) N(31) {17}[M(21,0)]{20} E(21,24) N(24) {25}[M(27,0)]{26} E(27,28) N(28) {26,19,15,7}[M(0,7pi/4)] E(0,27) E(0,29) N(29) {16}[M(20,0)]{17} E(20,21) N(21) {23}[M(26,0)]{25} E(26,27) N(27) {22}[M(25,0)]{23} E(25,26) N(26) {15}[M(17,7pi/4)]{16} E(17,20) N(20) {19}[M(23,pi/4)]{22} E(23,25) N(25) {14}[M(16,0)]{15} E(16,0) E(16,17) N(17) {13}[M(15,0)]{14} E(15,16) N(16) {18}[M(22,0)]{19} E(22,23) N(23) E(22,0) {12}[M(14,0)]{13} E(14,15) N(15) {1}[M(19,0)]{18} E(19,22) N(22) {11}[M(13,pi/4)]{12} E(13,14) N(14) [M(18,0)]{1} E(18,19) N(19) {10}[M(12,0)]{11} E(12,13) N(13) {11,3}[M(1,pi/4)] E(1,12) E(1,18) N(18) {9}[M(11,0)]{10} E(11,12) N(12) {8}[M(10,0)]{9} E(10,11) N(11) {7}[M(9,7pi/4)]{8} E(9,10) N(10) {6}[M(8,0)]{7} E(8,9) N(9) E(8,0) {5}[M(7,0)]{6} E(7,8) N(8) {4}[M(6,0)]{5} E(6,7) N(7) {3}[M(5,pi/4)]{4} E(5,6) N(6) {2}[M(4,0)]{3} E(4,5) N(5) E(4,1) [M(3,0)]{2} E(3,4) N(4) M(2,0) E(2,3) N(3)' >>> from graphix.command import CommandKind >>> pattern.to_ascii(target={CommandKind.M}) - '{17}[M(19)]{7} {10}[M(15,-pi/4)]{11} {17,12}[M(14)]{13} {17,9}[M(11)]{10} {0,1,5,10,13}[M(7,-pi/4)]{17} {1}[M(13,-7pi/4)]{12} {8}[M(10,-7pi/4)]{9} {17}[M(12)]{1} {6}[M(9)]{8} {8,3}[M(1,-pi/4)] {5}[M(8,-pi/4)]{6} {17,4}[M(6)]{5} [M(17)]{0} {3}[M(5,-7pi/4)]{4} M(0) {2}[M(4)]{3} [M(3)]{2} M(2)' + '{27}[M(31,0)]{28} [M(29,0)]{0} {26}[M(28,3pi/2)]{27} {17}[M(21,0)]{20} {25}[M(27,0)]{26} {26,19,15,7}[M(0,7pi/4)] {16}[M(20,0)]{17} {23}[M(26,0)]{25} {22}[M(25,0)]{23} {15}[M(17,7pi/4)]{16} {19}[M(23,pi/4)]{22} {14}[M(16,0)]{15} {13}[M(15,0)]{14} {18}[M(22,0)]{19} {12}[M(14,0)]{13} {1}[M(19,0)]{18} {11}[M(13,pi/4)]{12} [M(18,0)]{1} {10}[M(12,0)]{11} {11,3}[M(1,pi/4)] {9}[M(11,0)]{10} {8}[M(10,0)]{9} {7}[M(9,7pi/4)]{8} {6}[M(8,0)]{7} {5}[M(7,0)]{6} {4}[M(6,0)]{5} {3}[M(5,pi/4)]{4} {2}[M(4,0)]{3} [M(3,0)]{2} M(2,0)' >>> pattern.to_ascii(target={CommandKind.X, CommandKind.Z}) - 'Z(16,{11}) Z(18,{13}) Z(20,{7}) X(16,{15}) X(18,{14}) X(20,{19})' + 'X(24,{21}) Z(24,{20}) X(32,{31}) Z(32,{28}) X(30,{29}) Z(30,{0})' """ return pattern_to_str(self, OutputFormat.ASCII, left_to_right, limit, target) @@ -426,13 +426,14 @@ def to_latex( >>> circuit.ccx(0, 1, 2) >>> pattern = circuit.transpile().pattern >>> pattern.to_latex() - '\\(Z_{16}^{11}\\,Z_{18}^{13}\\,Z_{20}^{7}\\,X_{16}^{15}\\,X_{18}^{14}\\,X_{20}^{19}\\,{}_{17}[M_{19}]^{7}\\,{}_{10}[M_{15}^{-\\frac{\\pi}{4}}]^{11}\\,{}_{17,12}[M_{14}]^{13}\\,{}_{17,9}[M_{11}]^{10}\\,{}_{0,1,5,10,13}[M_{7}^{-\\frac{\\pi}{4}}]^{17}\\,{}_{1}[M_{13}^{-\\frac{7\\pi}{4}}]^{12}\\,{}_{8}[M_{10}^{-\\frac{7\\pi}{4}}]^{9}\\,{}_{17}[M_{12}]^{1}\\,{}_{6}[M_{9}]^{8}\\,{}_{8,3}[M_{1}^{-\\frac{\\pi}{4}}]\\,{}_{5}[M_{8}^{-\\frac{\\pi}{4}}]^{6}\\,{}_{17,4}[M_{6}]^{5}\\,[M_{17}]^{0}\\,{}_{3}[M_{5}^{-\\frac{7\\pi}{4}}]^{4}\\,M_{0}\\,{}_{2}[M_{4}]^{3}\\,[M_{3}]^{2}\\,M_{2}\\,E_{19,20}\\,E_{15,16}\\,E_{14,18}\\,E_{11,15}\\,E_{7,19}\\,E_{7,14}\\,E_{7,11}\\,E_{13,14}\\,E_{10,11}\\,E_{12,13}\\,E_{12,7}\\,E_{9,10}\\,E_{1,12}\\,E_{1,9}\\,E_{8,9}\\)...(27 more commands)' + '\\(X_{24}^{21}\\,Z_{24}^{20}\\,X_{32}^{31}\\,Z_{32}^{28}\\,X_{30}^{29}\\,Z_{30}^{0}\\,{}_{27}[M_{31}^{0}]^{28}\\,E_{31,32}\\,N_{32}\\,[M_{29}^{0}]^{0}\\,E_{29,30}\\,N_{30}\\,{}_{26}[M_{28}^{\\frac{3\\pi}{2}}]^{27}\\,E_{28,31}\\,N_{31}\\,{}_{17}[M_{21}^{0}]^{20}\\,E_{21,24}\\,N_{24}\\,{}_{25}[M_{27}^{0}]^{26}\\,E_{27,28}\\,N_{28}\\,{}_{26,19,15,7}[M_{0}^{\\frac{7\\pi}{4}}]\\,E_{0,27}\\,E_{0,29}\\,N_{29}\\,{}_{16}[M_{20}^{0}]^{17}\\,E_{20,21}\\,N_{21}\\,{}_{23}[M_{26}^{0}]^{25}\\,E_{26,27}\\,N_{27}\\,{}_{22}[M_{25}^{0}]^{23}\\,E_{25,26}\\,N_{26}\\,{}_{15}[M_{17}^{\\frac{7\\pi}{4}}]^{16}\\,E_{17,20}\\,N_{20}\\,{}_{19}[M_{23}^{\\frac{\\pi}{4}}]^{22}\\,E_{23,25}\\)...(63 more commands)' >>> pattern.to_latex(left_to_right=True) - '\\(N_{3}\\,N_{4}\\,N_{5}\\,N_{6}\\,N_{7}\\,N_{8}\\,N_{9}\\,N_{10}\\,N_{11}\\,N_{12}\\,N_{13}\\,N_{14}\\,N_{15}\\,N_{16}\\,N_{17}\\,N_{18}\\,N_{19}\\,N_{20}\\,E_{2,3}\\,E_{3,4}\\,E_{4,5}\\,E_{4,1}\\,E_{0,17}\\,E_{5,6}\\,E_{17,7}\\,E_{6,8}\\,E_{6,7}\\,E_{8,9}\\,E_{1,9}\\,E_{1,12}\\,E_{9,10}\\,E_{12,7}\\,E_{12,13}\\,E_{10,11}\\,E_{13,14}\\,E_{7,11}\\,E_{7,14}\\,E_{7,19}\\,E_{11,15}\\)...(27 more commands)' + '\\(N_{3}\\,E_{2,3}\\,M_{2}^{0}\\,N_{4}\\,E_{3,4}\\,[M_{3}^{0}]^{2}\\,E_{4,1}\\,N_{5}\\,E_{4,5}\\,{}_{2}[M_{4}^{0}]^{3}\\,N_{6}\\,E_{5,6}\\,{}_{3}[M_{5}^{\\frac{\\pi}{4}}]^{4}\\,N_{7}\\,E_{6,7}\\,{}_{4}[M_{6}^{0}]^{5}\\,N_{8}\\,E_{7,8}\\,{}_{5}[M_{7}^{0}]^{6}\\,E_{8,0}\\,N_{9}\\,E_{8,9}\\,{}_{6}[M_{8}^{0}]^{7}\\,N_{10}\\,E_{9,10}\\,{}_{7}[M_{9}^{\\frac{7\\pi}{4}}]^{8}\\,N_{11}\\,E_{10,11}\\,{}_{8}[M_{10}^{0}]^{9}\\,N_{12}\\,E_{11,12}\\,{}_{9}[M_{11}^{0}]^{10}\\,N_{18}\\,E_{1,18}\\,E_{1,12}\\,{}_{11,3}[M_{1}^{\\frac{\\pi}{4}}]\\,N_{13}\\,E_{12,13}\\,{}_{10}[M_{12}^{0}]^{11}\\)...(63 more commands)' + >>> from graphix.command import CommandKind >>> pattern.to_latex(target={CommandKind.M}) - '\\({}_{17}[M_{19}]^{7}\\,{}_{10}[M_{15}^{-\\frac{\\pi}{4}}]^{11}\\,{}_{17,12}[M_{14}]^{13}\\,{}_{17,9}[M_{11}]^{10}\\,{}_{0,1,5,10,13}[M_{7}^{-\\frac{\\pi}{4}}]^{17}\\,{}_{1}[M_{13}^{-\\frac{7\\pi}{4}}]^{12}\\,{}_{8}[M_{10}^{-\\frac{7\\pi}{4}}]^{9}\\,{}_{17}[M_{12}]^{1}\\,{}_{6}[M_{9}]^{8}\\,{}_{8,3}[M_{1}^{-\\frac{\\pi}{4}}]\\,{}_{5}[M_{8}^{-\\frac{\\pi}{4}}]^{6}\\,{}_{17,4}[M_{6}]^{5}\\,[M_{17}]^{0}\\,{}_{3}[M_{5}^{-\\frac{7\\pi}{4}}]^{4}\\,M_{0}\\,{}_{2}[M_{4}]^{3}\\,[M_{3}]^{2}\\,M_{2}\\)' + '\\({}_{27}[M_{31}^{0}]^{28}\\,[M_{29}^{0}]^{0}\\,{}_{26}[M_{28}^{\\frac{3\\pi}{2}}]^{27}\\,{}_{17}[M_{21}^{0}]^{20}\\,{}_{25}[M_{27}^{0}]^{26}\\,{}_{26,19,15,7}[M_{0}^{\\frac{7\\pi}{4}}]\\,{}_{16}[M_{20}^{0}]^{17}\\,{}_{23}[M_{26}^{0}]^{25}\\,{}_{22}[M_{25}^{0}]^{23}\\,{}_{15}[M_{17}^{\\frac{7\\pi}{4}}]^{16}\\,{}_{19}[M_{23}^{\\frac{\\pi}{4}}]^{22}\\,{}_{14}[M_{16}^{0}]^{15}\\,{}_{13}[M_{15}^{0}]^{14}\\,{}_{18}[M_{22}^{0}]^{19}\\,{}_{12}[M_{14}^{0}]^{13}\\,{}_{1}[M_{19}^{0}]^{18}\\,{}_{11}[M_{13}^{\\frac{\\pi}{4}}]^{12}\\,[M_{18}^{0}]^{1}\\,{}_{10}[M_{12}^{0}]^{11}\\,{}_{11,3}[M_{1}^{\\frac{\\pi}{4}}]\\,{}_{9}[M_{11}^{0}]^{10}\\,{}_{8}[M_{10}^{0}]^{9}\\,{}_{7}[M_{9}^{\\frac{7\\pi}{4}}]^{8}\\,{}_{6}[M_{8}^{0}]^{7}\\,{}_{5}[M_{7}^{0}]^{6}\\,{}_{4}[M_{6}^{0}]^{5}\\,{}_{3}[M_{5}^{\\frac{\\pi}{4}}]^{4}\\,{}_{2}[M_{4}^{0}]^{3}\\,[M_{3}^{0}]^{2}\\,M_{2}^{0}\\)' >>> pattern.to_latex(target={CommandKind.X, CommandKind.Z}) - '\\(Z_{16}^{11}\\,Z_{18}^{13}\\,Z_{20}^{7}\\,X_{16}^{15}\\,X_{18}^{14}\\,X_{20}^{19}\\)' + '\\(X_{24}^{21}\\,Z_{24}^{20}\\,X_{32}^{31}\\,Z_{32}^{28}\\,X_{30}^{29}\\,Z_{30}^{0}\\)' """ return pattern_to_str(self, OutputFormat.LaTeX, left_to_right, limit, target) @@ -461,13 +462,14 @@ def to_unicode( >>> circuit.ccx(0, 1, 2) >>> pattern = circuit.transpile().pattern >>> pattern.to_unicode() - 'Z₁₆¹¹ Z₁₈¹³ Z₂₀⁷ X₁₆¹⁵ X₁₈¹⁴ X₂₀¹⁹ ₁₇[M₁₉]⁷ ₁₀[M₁₅(-π/4)]¹¹ ₁₇₊₁₂[M₁₄]¹³ ₁₇₊₉[M₁₁]¹⁰ ₀₊₁₊₅₊₁₀₊₁₃[M₇(-π/4)]¹⁷ ₁[M₁₃(-7π/4)]¹² ₈[M₁₀(-7π/4)]⁹ ₁₇[M₁₂]¹ ₆[M₉]⁸ ₈₊₃[M₁(-π/4)] ₅[M₈(-π/4)]⁶ ₁₇₊₄[M₆]⁵ [M₁₇]⁰ ₃[M₅(-7π/4)]⁴ M₀ ₂[M₄]³ [M₃]² M₂ E₁₉₋₂₀ E₁₅₋₁₆ E₁₄₋₁₈ E₁₁₋₁₅ E₇₋₁₉ E₇₋₁₄ E₇₋₁₁ E₁₃₋₁₄ E₁₀₋₁₁ E₁₂₋₁₃ E₁₂₋₇ E₉₋₁₀ E₁₋₁₂ E₁₋₉ E₈₋₉...(27 more commands)' + 'X₂₄²¹ Z₂₄²⁰ X₃₂³¹ Z₃₂²⁸ X₃₀²⁹ Z₃₀⁰ ₂₇[M₃₁(0)]²⁸ E₃₁₋₃₂ N₃₂ [M₂₉(0)]⁰ E₂₉₋₃₀ N₃₀ ₂₆[M₂₈(3π/2)]²⁷ E₂₈₋₃₁ N₃₁ ₁₇[M₂₁(0)]²⁰ E₂₁₋₂₄ N₂₄ ₂₅[M₂₇(0)]²⁶ E₂₇₋₂₈ N₂₈ ₂₆₊₁₉₊₁₅₊₇[M₀(7π/4)] E₀₋₂₇ E₀₋₂₉ N₂₉ ₁₆[M₂₀(0)]¹⁷ E₂₀₋₂₁ N₂₁ ₂₃[M₂₆(0)]²⁵ E₂₆₋₂₇ N₂₇ ₂₂[M₂₅(0)]²³ E₂₅₋₂₆ N₂₆ ₁₅[M₁₇(7π/4)]¹⁶ E₁₇₋₂₀ N₂₀ ₁₉[M₂₃(π/4)]²² E₂₃₋₂₅...(63 more commands)' >>> pattern.to_unicode(left_to_right=True) - 'N₃ N₄ N₅ N₆ N₇ N₈ N₉ N₁₀ N₁₁ N₁₂ N₁₃ N₁₄ N₁₅ N₁₆ N₁₇ N₁₈ N₁₉ N₂₀ E₂₋₃ E₃₋₄ E₄₋₅ E₄₋₁ E₀₋₁₇ E₅₋₆ E₁₇₋₇ E₆₋₈ E₆₋₇ E₈₋₉ E₁₋₉ E₁₋₁₂ E₉₋₁₀ E₁₂₋₇ E₁₂₋₁₃ E₁₀₋₁₁ E₁₃₋₁₄ E₇₋₁₁ E₇₋₁₄ E₇₋₁₉ E₁₁₋₁₅...(27 more commands)' + 'N₃ E₂₋₃ M₂(0) N₄ E₃₋₄ [M₃(0)]² E₄₋₁ N₅ E₄₋₅ ₂[M₄(0)]³ N₆ E₅₋₆ ₃[M₅(π/4)]⁴ N₇ E₆₋₇ ₄[M₆(0)]⁵ N₈ E₇₋₈ ₅[M₇(0)]⁶ E₈₋₀ N₉ E₈₋₉ ₆[M₈(0)]⁷ N₁₀ E₉₋₁₀ ₇[M₉(7π/4)]⁸ N₁₁ E₁₀₋₁₁ ₈[M₁₀(0)]⁹ N₁₂ E₁₁₋₁₂ ₉[M₁₁(0)]¹⁰ N₁₈ E₁₋₁₈ E₁₋₁₂ ₁₁₊₃[M₁(π/4)] N₁₃ E₁₂₋₁₃ ₁₀[M₁₂(0)]¹¹...(63 more commands)' + >>> from graphix.command import CommandKind >>> pattern.to_unicode(target={CommandKind.M}) - '₁₇[M₁₉]⁷ ₁₀[M₁₅(-π/4)]¹¹ ₁₇₊₁₂[M₁₄]¹³ ₁₇₊₉[M₁₁]¹⁰ ₀₊₁₊₅₊₁₀₊₁₃[M₇(-π/4)]¹⁷ ₁[M₁₃(-7π/4)]¹² ₈[M₁₀(-7π/4)]⁹ ₁₇[M₁₂]¹ ₆[M₉]⁸ ₈₊₃[M₁(-π/4)] ₅[M₈(-π/4)]⁶ ₁₇₊₄[M₆]⁵ [M₁₇]⁰ ₃[M₅(-7π/4)]⁴ M₀ ₂[M₄]³ [M₃]² M₂' + '₂₇[M₃₁(0)]²⁸ [M₂₉(0)]⁰ ₂₆[M₂₈(3π/2)]²⁷ ₁₇[M₂₁(0)]²⁰ ₂₅[M₂₇(0)]²⁶ ₂₆₊₁₉₊₁₅₊₇[M₀(7π/4)] ₁₆[M₂₀(0)]¹⁷ ₂₃[M₂₆(0)]²⁵ ₂₂[M₂₅(0)]²³ ₁₅[M₁₇(7π/4)]¹⁶ ₁₉[M₂₃(π/4)]²² ₁₄[M₁₆(0)]¹⁵ ₁₃[M₁₅(0)]¹⁴ ₁₈[M₂₂(0)]¹⁹ ₁₂[M₁₄(0)]¹³ ₁[M₁₉(0)]¹⁸ ₁₁[M₁₃(π/4)]¹² [M₁₈(0)]¹ ₁₀[M₁₂(0)]¹¹ ₁₁₊₃[M₁(π/4)] ₉[M₁₁(0)]¹⁰ ₈[M₁₀(0)]⁹ ₇[M₉(7π/4)]⁸ ₆[M₈(0)]⁷ ₅[M₇(0)]⁶ ₄[M₆(0)]⁵ ₃[M₅(π/4)]⁴ ₂[M₄(0)]³ [M₃(0)]² M₂(0)' >>> pattern.to_unicode(target={CommandKind.X, CommandKind.Z}) - 'Z₁₆¹¹ Z₁₈¹³ Z₂₀⁷ X₁₆¹⁵ X₁₈¹⁴ X₂₀¹⁹' + 'X₂₄²¹ Z₂₄²⁰ X₃₂³¹ Z₃₂²⁸ X₃₀²⁹ Z₃₀⁰' """ return pattern_to_str(self, OutputFormat.Unicode, left_to_right, limit, target) @@ -531,7 +533,7 @@ def shift_signals(self, method: str = "direct") -> dict[int, set[int]]: ---------- method : str, optional 'direct' shift_signals is executed on a conventional Pattern sequence. - 'mc' shift_signals is done using the original algorithm on the measurement calculus paper. + 'mc' shift_signals is done using the original algorithm in the measurement calculus paper. Returns ------- diff --git a/graphix/qasm3_exporter.py b/graphix/qasm3_exporter.py index 7e4a9336d..19d152f64 100644 --- a/graphix/qasm3_exporter.py +++ b/graphix/qasm3_exporter.py @@ -22,25 +22,47 @@ from graphix.instruction import InstructionType -def circuit_to_qasm3(circuit: Circuit) -> str: - """Export circuit instructions to OpenQASM 3.0 representation. +def circuit_to_qasm3(circuit: Circuit, *, transpile: bool = True) -> str: + """Export circuit instructions to an OpenQASM 3.0 string. + + Parameters + ---------- + circuit: Circuit + The quantum circuit to be exported. + + transpile: bool, optional + If ``True`` (the default), the circuit is first transpiled so that any J gates + are replaced by equivalent RZ and H gates. + If ``False``, a ``ValueError`` is raised when the circuit contains a J gate. Returns ------- str - The OpenQASM 3.0 string representation of the circuit. + The OpenQASM 3.0 representation of the supplied circuit. """ - return "\n".join(circuit_to_qasm3_lines(circuit)) + return "\n".join(circuit_to_qasm3_lines(circuit, transpile=transpile)) -def circuit_to_qasm3_lines(circuit: Circuit) -> Iterator[str]: +def circuit_to_qasm3_lines(circuit: Circuit, *, transpile: bool = True) -> Iterator[str]: """Export circuit instructions to line-by-line OpenQASM 3.0 representation. + Parameters + ---------- + circuit: Circuit + The quantum circuit to be exported. + + transpile: bool, optional + If ``True`` (the default), the circuit is first transpiled so that any J gates + are replaced by equivalent RZ and H gates. + If ``False``, a ``ValueError`` is raised when the circuit contains a J gate. + Returns ------- Iterator[str] - The OpenQASM 3.0 lines that represent the circuit. + An iterator over the OpenQASM 3.0 lines that represent the circuit. """ + if transpile: + circuit = circuit.transpile_j_to_rzh().transpile_measurements_to_z_axis() yield "OPENQASM 3;" yield 'include "stdgates.inc";' yield f"qubit[{circuit.width}] q;" @@ -72,12 +94,28 @@ def angle_to_qasm3(angle: ParameterizedAngle) -> str: def instruction_to_qasm3(instruction: InstructionType) -> str: - """Get the OpenQASM3 representation of a single circuit instruction.""" + """Get the OpenQASM3 representation of a single circuit instruction. + + Parameters + ---------- + instruction : Instruction + The instruction to convert. + + Returns + ------- + str + The OpenQASM3 representation of the instruction. + + Raises + ------ + ValueError + If the instruction is not supported by OpenQASM3. + """ match instruction.kind: case InstructionKind.M: if instruction.axis != Axis.Z: raise ValueError( - "OpenQASM3 only supports measurements on Z axis. Use `Circuit.transpile_measurements_to_z_axis` to rewrite measurements on X and Y axes." + "OpenQASM3 only supports measurements on Z axis. Use `Circuit.transpile_measurements_to_z_axis` to rewrite measurements on X and Y axes, or setting `transpile=True`." ) return f"b[{instruction.target}] = measure q[{instruction.target}]" case InstructionKind.RX | InstructionKind.RY | InstructionKind.RZ: @@ -85,6 +123,10 @@ def instruction_to_qasm3(instruction: InstructionType) -> str: return qasm3_gate_call( instruction.kind.name.lower(), args=[angle], operands=[qasm3_qubit(instruction.target)] ) + case InstructionKind.J: + raise ValueError( + "J gates must be decomposed before QASM3 export using `Circuit.transpile_j_to_rzh`, or setting `transpile=True`." + ) case InstructionKind.H | InstructionKind.S | InstructionKind.X | InstructionKind.Y | InstructionKind.Z: return qasm3_gate_call(instruction.kind.name.lower(), [qasm3_qubit(instruction.target)]) case InstructionKind.I: diff --git a/graphix/random_objects.py b/graphix/random_objects.py index d5c89f788..c966c4ae3 100644 --- a/graphix/random_objects.py +++ b/graphix/random_objects.py @@ -377,6 +377,7 @@ def rand_circuit( depth: int, rng: Generator | None = None, *, + use_j: bool = True, use_cz: bool = True, use_rzz: bool = False, use_ccx: bool = False, @@ -393,6 +394,8 @@ def rand_circuit( Number of alternating entangling and single-qubit layers. rng : numpy.random.Generator, optional Random number generator. A default generator is created if ``None``. + use_j : bool, optional + If ``True`` add J gates in each layer (default: ``False``). use_cz : bool, optional If ``True`` add CZ gates in each layer (default: ``True``). use_rzz : bool, optional @@ -424,6 +427,7 @@ def rand_circuit( circuit.z, circuit.y, *parametric_gate_choice, + *((functools.partial(circuit.j, angle=ANGLE_PI / 4),) if use_j else ()), ) for _ in range(depth): for j, k in _genpair(nqubits, 2, rng): @@ -437,7 +441,7 @@ def rand_circuit( if use_ccx: for j, k, l in _gentriplet(nqubits, 2, rng): circuit.ccx(j, k, l) - for j, k in _genpair(nqubits, 4, rng): + for j, k in _genpair(nqubits, 2, rng): circuit.swap(j, k) for j in range(nqubits): ind = rng.integers(len(gate_choice)) diff --git a/graphix/sim/base_backend.py b/graphix/sim/base_backend.py index ef1977b3c..5dea50ce1 100644 --- a/graphix/sim/base_backend.py +++ b/graphix/sim/base_backend.py @@ -230,8 +230,7 @@ class NodeIndex: def __init__(self) -> None: """Initialize an empty mapping between nodes and qubit indices.""" - self.__dict = {} - self.__list = [] + self.clear() def __getitem__(self, index: int) -> int: """Return the qubit node associated with the specified index. @@ -316,6 +315,11 @@ def swap(self, i: int, j: int) -> None: self.__dict[node_i] = j self.__dict[node_j] = i + def clear(self) -> None: + """Delete all the elements, restoring the index to its initial state.""" + self.__dict = {} + self.__list = [] + class NoiseNotSupportedError(Exception): """Exception raised when `apply_channel` is called on a backend that does not support noise.""" @@ -443,6 +447,17 @@ def swap(self, qubits: tuple[int, int]) -> None: (control, target) qubit indices """ + @abstractmethod + def permute(self, permutation: Sequence[int]) -> None: + """Reorder the qubits. + + Parameters + ---------- + permutation: Sequence[int] + The permutation to apply. For each position in the resulting order, + the value gives the index of the qubit in the original ordering. + """ + def apply_noise(self, qubits: Sequence[int], noise: Noise) -> None: # noqa: ARG002 """Apply noise. @@ -577,6 +592,11 @@ class Backend(Generic[_StateT_co]): # (specifically, as a parameter of `__init__`) since `_StateT_co` is covariant. state: _StateT_co = dataclasses.field(init=False) + @property + @abstractmethod + def nqubit(self) -> int: + """Return the number of qubits.""" + @abstractmethod def add_nodes(self, nodes: Sequence[int], data: Data = BasicStates.PLUS) -> None: r""" @@ -659,7 +679,7 @@ def entangle_nodes(self, edge: tuple[int, int]) -> None: """ @abstractmethod - def finalize(self, output_nodes: Iterable[int]) -> None: + def finalize(self, output_nodes: Sequence[int]) -> None: """To be run at the end of pattern simulation to convey the order of output nodes.""" @abstractmethod @@ -818,6 +838,13 @@ def apply_clifford(self, node: int, clifford: Clifford) -> None: loc = self.node_index.index(node) self.state.evolve_single(clifford.matrix, loc) + @override + def finalize(self, output_nodes: Sequence[int]) -> None: + """To be run at the end of pattern simulation.""" + self.state.permute([self.node_index.index(node) for node in output_nodes]) + self.node_index.clear() + self.node_index.extend(output_nodes) + def sort_qubits(self, output_nodes: Iterable[int]) -> None: """Sort the qubit order in internal statevector.""" for i, ind in enumerate(output_nodes): @@ -826,12 +853,18 @@ def sort_qubits(self, output_nodes: Iterable[int]) -> None: self.state.swap((i, move_from)) self.node_index.swap(i, move_from) - @override - def finalize(self, output_nodes: Iterable[int]) -> None: - """To be run at the end of pattern simulation.""" - self.sort_qubits(output_nodes) + # @override + # def finalize(self, output_nodes: Iterable[int]) -> None: + # """To be run at the end of pattern simulation.""" + # from copy import copy + # duplicate = copy(self) + # duplicate.sort_qubits(output_nodes) + # self.state.permute([self.node_index.index(node) for node in output_nodes]) + # if not self.state.isclose(duplicate.state): + # assert False, output_nodes @property + @override def nqubit(self) -> int: """Return the number of qubits of the current state.""" return self.state.nqubit diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index cb5570e2a..29538cc79 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -20,7 +20,7 @@ from graphix.channels import KrausChannel from graphix.parameter import Expression, ExpressionOrFloat, ExpressionOrSupportsComplex from graphix.sim.base_backend import DenseState, DenseStateBackend, Matrix, kron, matmul, outer, tensordot, vdot -from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, Statevec +from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, Statevec, _check_permutation from graphix.states import BasicStates, State if TYPE_CHECKING: @@ -137,6 +137,7 @@ def add_nodes(self, nqubit: int, data: Data) -> None: - A multi-qubit state vector of dimension :math:`2^n` initializes the new nodes jointly. - A density matrix must have shape :math:`2^n \times 2^n`, and is used to jointly initialize the new nodes. + - The type of nodes to be added is inferred from the type of the existing ``DensityMatrix``. Notes ----- @@ -259,6 +260,8 @@ def tensor(self, other: DensityMatrix) -> None: """ if not isinstance(other, DensityMatrix): other = DensityMatrix(other) + if self.rho.dtype == np.object_ and other.rho.dtype != np.object_: + other.rho = other.rho.astype(np.object_, copy=False) # pragma: nocover self.rho = kron(self.rho, other.rho) def cnot(self, edge: tuple[int, int]) -> None: @@ -282,6 +285,17 @@ def swap(self, qubits: tuple[int, int]) -> None: """ self.evolve(SWAP_TENSOR.reshape(4, 4), qubits) + @override + def permute(self, permutation: Sequence[int]) -> None: + nqubit = self.nqubit + _check_permutation(permutation, nqubit) + tensor_shape = [2] * (2 * nqubit) + perm_cols = [i + nqubit for i in permutation] + full_permutation = [*permutation, *perm_cols] + rho_tensor = self.rho.reshape(tensor_shape) + rho_permuted_tensor = np.transpose(rho_tensor, axes=full_permutation) + self.rho = rho_permuted_tensor.reshape((2**nqubit, 2**nqubit)) + def entangle(self, edge: tuple[int, int]) -> None: """Connect graph nodes. diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 1f7dd23e6..2ee79149a 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -165,6 +165,7 @@ def add_nodes(self, nqubit: int, data: Data) -> None: - A single-qubit state vector will be broadcast to all nodes. - A multi-qubit state vector of dimension :math:`2^n`, where :math:`n = \mathrm{len}(nodes)`, initializes the new nodes jointly. + - The type of nodes to be added is inferred from the type of the existing ``Statevec``. Notes ----- @@ -301,6 +302,8 @@ def tensor(self, other: Statevec) -> None: psi_other = other.psi.flatten() total_num = len(self.dims()) + len(other.dims()) + if psi_self.dtype == np.object_ and psi_other.dtype != np.object_: + psi_other = psi_other.astype(np.object_, copy=False) # pragma: nocover self.psi = kron(psi_self, psi_other).reshape((2,) * total_num) def cnot(self, qubits: tuple[int, int]) -> None: @@ -330,6 +333,11 @@ def swap(self, qubits: tuple[int, int]) -> None: # sort back axes self.psi = np.moveaxis(psi, (0, 1), qubits) + @override + def permute(self, permutation: Sequence[int]) -> None: + _check_permutation(permutation, self.nqubit) + self.psi = np.transpose(self.psi, permutation) + def normalize(self) -> None: """Normalize the state in-place.""" # Note that the following calls to `astype` are guaranteed to @@ -582,3 +590,10 @@ def _format_encoding(nqubit: int, i: int, encoding: _ENCODING) -> str: if encoding == "LSB": return output[::-1] return output + + +def _check_permutation(permutation: Sequence[int], nqubits: int) -> None: + if len(permutation) != nqubits: + raise ValueError(f"Permutation has length {len(permutation)}, but {nqubits} qubits expected.") + if set(permutation) != set(range(nqubits)): + raise ValueError(f"{permutation} is not a permutation.") diff --git a/graphix/sim/tensornet.py b/graphix/sim/tensornet.py index f9e8c56f7..f438bd0c6 100644 --- a/graphix/sim/tensornet.py +++ b/graphix/sim/tensornet.py @@ -801,6 +801,12 @@ def apply_clifford(self, node: int, clifford: Clifford) -> None: def finalize(self, output_nodes: Iterable[int]) -> None: """Do nothing.""" + @property + @override + def nqubit(self) -> int: + """Raise NotImplementedError: the current number of qubits is not well-defined in the current implementation of the TensorNetworkBackend.""" + raise NotImplementedError + def gen_str() -> str: """Generate dummy string for einsum.""" diff --git a/graphix/simulator.py b/graphix/simulator.py index e81439ebc..50410d1e2 100644 --- a/graphix/simulator.py +++ b/graphix/simulator.py @@ -407,6 +407,25 @@ def run( Stack level to use for warnings. Defaults to 1, meaning that warnings are reported at this function's call site. """ + # Check whether the backend is properly initialized. We + # disable the check for TensorNetworkBackend because its + # current behavior differs from other backends. + if not isinstance(self.backend, TensorNetworkBackend): + initial_nqubit = self.backend.nqubit + if input_state is None: + # No explicit state supplied: the backend must already contain the + # required input qubits. + input_nodes_len = len(self.pattern.input_nodes) + if initial_nqubit != input_nodes_len: + raise ValueError( + f"`input_state` is `None`: the backend is expected to have {input_nodes_len} input nodes already prepared, but {initial_nqubit} were found." + ) + # An explicit state was supplied: the backend must start with a clean + # state (no pre-allocated qubits). + elif initial_nqubit != 0: + raise ValueError( + f"`input_state` is not `None`: the backend is expected to have no pre-allocated qubits, but has {initial_nqubit} qubits." + ) if input_state is not None: self.backend.add_nodes(self.pattern.input_nodes, input_state) if self.noise_model is None: diff --git a/graphix/transpiler.py b/graphix/transpiler.py index 7859d3f0f..aa49d76bb 100644 --- a/graphix/transpiler.py +++ b/graphix/transpiler.py @@ -6,24 +6,30 @@ from __future__ import annotations +import enum from dataclasses import dataclass +from enum import Enum from typing import TYPE_CHECKING, Generic, SupportsFloat, TypeVar, overload +import networkx as nx + # assert_never introduced in Python 3.11 # override introduced in Python 3.12 from typing_extensions import assert_never, override from graphix import command, instruction, parameter from graphix.branch_selector import BranchSelector, RandomBranchSelector -from graphix.command import E, M, N, X, Z +from graphix.flow.core import CausalFlow, _corrections_to_partial_order_layers from graphix.fundamentals import ANGLE_PI, Axis from graphix.instruction import InstructionKind, InstructionVisitor -from graphix.measurements import Measurement, PauliMeasurement +from graphix.measurements import BlochMeasurement, Measurement, Outcome, PauliMeasurement +from graphix.opengraph import OpenGraph from graphix.ops import Ops +from graphix.optimization import StandardizedPattern from graphix.pattern import Pattern from graphix.sim.base_backend import DenseStateBackend from graphix.sim.density_matrix import DensityMatrixBackend -from graphix.sim.statevec import StatevectorBackend +from graphix.sim.statevec import Statevec, StatevectorBackend if TYPE_CHECKING: from collections.abc import Callable, Iterable, Iterator, Mapping, Sequence @@ -31,33 +37,45 @@ from numpy.random import Generator - from graphix.command import CommandType + from graphix.command import Node from graphix.fundamentals import ParameterizedAngle - from graphix.instruction import InstructionType, InstructionTypeWithoutRZZ + from graphix.instruction import InstructionType from graphix.parameter import ExpressionOrFloat, Parameter + from graphix.pattern import Pattern from graphix.sim import Data from graphix.sim.base_backend import DenseState, Matrix from graphix.sim.density_matrix import DensityMatrix - from graphix.sim.statevec import Statevec _BuiltinDenseStateBackend = DensityMatrixBackend | StatevectorBackend _DenseStateBackendLiteral = Literal["statevector", "densitymatrix"] -_DenseStateT_co = TypeVar("_DenseStateT_co", bound="DenseState", covariant=True) _DenseStateT = TypeVar("_DenseStateT", bound="DenseState") -@dataclass(frozen=True) -class TranspileResult: - """ - The result of a transpilation. - - pattern : :class:`graphix.pattern.Pattern` object - classical_outputs : tuple[int,...], index of nodes measured with *M* gates - """ +@dataclass(frozen=True, slots=True) +class TranspiledPattern: + """A transpiled pattern.""" pattern: Pattern - classical_outputs: tuple[int, ...] + + classical_outputs: tuple[Node, ...] + """Nodes measured with circuit measurements, in the order of the gates.""" + + +@dataclass(frozen=True, slots=True) +class TranspiledFlow: + """A transpiled causal flow.""" + + flow: CausalFlow[BlochMeasurement] + + classical_outputs: dict[int, command.M] + """M commands for nodes measured with circuit measurements.""" + + def to_pattern(self) -> TranspiledPattern: + """Return the transpiled pattern.""" + pattern = StandardizedPattern.from_pattern(self.flow.to_corrections().to_pattern()).to_space_optimal_pattern() + pattern.extend(self.classical_outputs.values()) + return TranspiledPattern(pattern, tuple(self.classical_outputs.keys())) @dataclass(frozen=True) @@ -157,6 +175,8 @@ def add(self, instr: InstructionType) -> None: self.ry(instr.target, instr.angle) case InstructionKind.RZ: self.rz(instr.target, instr.angle) + case InstructionKind.J: + self.j(instr.target, instr.angle) case _: assert_never(instr.kind) @@ -308,6 +328,19 @@ def rz(self, qubit: int, angle: ParameterizedAngle) -> None: assert qubit in self.active_qubits self.instruction.append(instruction.RZ(target=qubit, angle=angle)) + def j(self, qubit: int, angle: ParameterizedAngle) -> None: + """Apply a J rotation gate. + + Parameters + ---------- + qubit : int + target qubit + angle : ParameterizedAngle + rotation angle in units of π + """ + assert qubit in self.active_qubits + self.instruction.append(instruction.J(target=qubit, angle=angle)) + def r(self, qubit: int, axis: Axis, angle: ParameterizedAngle) -> None: """Apply a rotation gate on the given axis. @@ -401,559 +434,95 @@ def m(self, qubit: int, axis: Axis) -> None: self.instruction.append(instruction.M(target=qubit, axis=axis)) self.active_qubits.remove(qubit) - def transpile(self) -> TranspileResult: - """Transpile the circuit to a pattern. - - Returns - ------- - result : :class:`TranspileResult` object - """ - n_node = self.width - out: list[int | None] = list(range(self.width)) - pattern = Pattern(input_nodes=list(range(self.width))) - classical_outputs = [] - for instr in _transpile_rzz(self.instruction): - match instr.kind: - case instruction.InstructionKind.CZ: - target0 = _check_target(out, instr.targets[0]) - target1 = _check_target(out, instr.targets[1]) - seq = self._cz_command(target0, target1) - pattern.extend(seq) - case instruction.InstructionKind.CNOT: - ancilla = [n_node, n_node + 1] - control = _check_target(out, instr.control) - target = _check_target(out, instr.target) - out[instr.control], out[instr.target], seq = self._cnot_command(control, target, ancilla) - pattern.extend(seq) - n_node += 2 - case instruction.InstructionKind.SWAP: - target0 = _check_target(out, instr.targets[0]) - target1 = _check_target(out, instr.targets[1]) - out[instr.targets[0]], out[instr.targets[1]] = ( - target1, - target0, - ) - case instruction.InstructionKind.I: - pass - case instruction.InstructionKind.H: - single_ancilla = n_node - target = _check_target(out, instr.target) - out[instr.target], seq = self._h_command(target, single_ancilla) - pattern.extend(seq) - n_node += 1 - case instruction.InstructionKind.S: - ancilla = [n_node, n_node + 1] - target = _check_target(out, instr.target) - out[instr.target], seq = self._s_command(target, ancilla) - pattern.extend(seq) - n_node += 2 - case instruction.InstructionKind.X: - ancilla = [n_node, n_node + 1] - target = _check_target(out, instr.target) - out[instr.target], seq = self._x_command(target, ancilla) - pattern.extend(seq) - n_node += 2 - case instruction.InstructionKind.Y: - ancilla = [n_node, n_node + 1, n_node + 2, n_node + 3] - target = _check_target(out, instr.target) - out[instr.target], seq = self._y_command(target, ancilla) - pattern.extend(seq) - n_node += 4 - case instruction.InstructionKind.Z: - ancilla = [n_node, n_node + 1] - target = _check_target(out, instr.target) - out[instr.target], seq = self._z_command(target, ancilla) - pattern.extend(seq) - n_node += 2 - case instruction.InstructionKind.RX: - ancilla = [n_node, n_node + 1] - target = _check_target(out, instr.target) - out[instr.target], seq = self._rx_command(target, ancilla, instr.angle) - pattern.extend(seq) - n_node += 2 - case instruction.InstructionKind.RY: - ancilla = [n_node, n_node + 1, n_node + 2, n_node + 3] - target = _check_target(out, instr.target) - out[instr.target], seq = self._ry_command(target, ancilla, instr.angle) - pattern.extend(seq) - n_node += 4 - case instruction.InstructionKind.RZ: - ancilla = [n_node, n_node + 1] - target = _check_target(out, instr.target) - out[instr.target], seq = self._rz_command(target, ancilla, instr.angle) - pattern.extend(seq) - n_node += 2 - case instruction.InstructionKind.CCX: - ancilla = [n_node + i for i in range(18)] - control0 = _check_target(out, instr.controls[0]) - control1 = _check_target(out, instr.controls[1]) - target = _check_target(out, instr.target) - ( - out[instr.controls[0]], - out[instr.controls[1]], - out[instr.target], - seq, - ) = self._ccx_command( - control0, - control1, - target, - ancilla, - ) - pattern.extend(seq) - n_node += 18 - case instruction.InstructionKind.M: - target = _check_target(out, instr.target) - seq = self._m_command(target, instr.axis) - pattern.extend(seq) - classical_outputs.append(target) - out[instr.target] = None - case _: - assert_never(instr.kind) - output_nodes = [node for node in out if node is not None] - pattern.reorder_output_nodes(output_nodes) - return TranspileResult(pattern, tuple(classical_outputs)) - - @classmethod - def _cnot_command( - cls, control_node: int, target_node: int, ancilla: Sequence[int] - ) -> tuple[int, int, list[command.CommandType]]: - """MBQC commands for CNOT gate. + def transpile_to_causal_flow(self) -> TranspiledFlow: + """Transpile a circuit via J-∧z decomposition to a causal flow. Parameters ---------- - control_node : int - control node on graph - target : int - target node on graph - ancilla : list of two ints - ancilla node indices to be added to graph + self: the circuit to transpile. Returns ------- - control_out : int - control node on graph after the gate - target_out : int - target node on graph after the gate - commands : list - list of MBQC commands + the result of the transpilation: a causal flow and classical outputs. """ - assert len(ancilla) == 2 - seq: list[CommandType] = [N(node=ancilla[0]), N(node=ancilla[1])] - seq.extend( - ( - E(nodes=(target_node, ancilla[0])), - E(nodes=(control_node, ancilla[0])), - E(nodes=(ancilla[0], ancilla[1])), - M(node=target_node), - M(node=ancilla[0], s_domain={target_node}), - X(node=ancilla[1], domain={ancilla[0]}), - Z(node=ancilla[1], domain={target_node}), - Z(node=control_node, domain={target_node}), - ) - ) - return control_node, ancilla[1], seq - - @classmethod - def _cz_command(cls, target_1: int, target_2: int) -> list[CommandType]: - """MBQC commands for CZ gate. - - Parameters - ---------- - target_1 : int - target node on graph - target_2 : int - other target node on graph - - Returns - ------- - commands : list - list of MBQC commands - """ - return [E(nodes=(target_1, target_2))] - - @classmethod - def _m_command(cls, input_node: int, axis: Axis) -> list[CommandType]: - """MBQC commands for measuring qubit. - - Parameters - ---------- - input_node : int - target node on graph - axis : Axis - measurement basis - - Returns - ------- - commands : list - list of MBQC commands - """ - # `measurement.angle` and `M.angle` are both expressed in units of π. - return [M(input_node, PauliMeasurement(axis))] - - @classmethod - def _h_command(cls, input_node: int, ancilla: int) -> tuple[int, list[CommandType]]: - """MBQC commands for Hadamard gate. - - Parameters - ---------- - input_node : int - target node on graph - ancilla : int - ancilla node index to be added - - Returns - ------- - out_node : int - control node on graph after the gate - commands : list - list of MBQC commands - """ - seq: list[CommandType] = [N(node=ancilla)] - seq.extend((E(nodes=(input_node, ancilla)), M(node=input_node), X(node=ancilla, domain={input_node}))) - return ancilla, seq - - @classmethod - def _s_command(cls, input_node: int, ancilla: Sequence[int]) -> tuple[int, list[command.CommandType]]: - """MBQC commands for S gate. - - Parameters - ---------- - input_node : int - input node index - ancilla : list of two ints - ancilla node indices to be added to graph - - Returns - ------- - out_node : int - control node on graph after the gate - commands : list - list of MBQC commands - """ - assert len(ancilla) == 2 - seq: list[CommandType] = [N(node=ancilla[0]), command.N(node=ancilla[1])] - seq.extend( - ( - E(nodes=(input_node, ancilla[0])), - E(nodes=(ancilla[0], ancilla[1])), - M(input_node, -Measurement.Y), - M(node=ancilla[0], s_domain={input_node}), - X(node=ancilla[1], domain={ancilla[0]}), - Z(node=ancilla[1], domain={input_node}), - ) - ) - return ancilla[1], seq - - @classmethod - def _x_command(cls, input_node: int, ancilla: Sequence[int]) -> tuple[int, list[command.CommandType]]: - """MBQC commands for Pauli X gate. - - Parameters - ---------- - input_node : int - input node index - ancilla : list of two ints - ancilla node indices to be added to graph - - Returns - ------- - out_node : int - control node on graph after the gate - commands : list - list of MBQC commands - """ - assert len(ancilla) == 2 - seq: list[CommandType] = [N(node=ancilla[0]), N(node=ancilla[1])] - seq.extend( - ( - E(nodes=(input_node, ancilla[0])), - E(nodes=(ancilla[0], ancilla[1])), - M(node=input_node), - M(ancilla[0], -Measurement.X, s_domain={input_node}), - X(node=ancilla[1], domain={ancilla[0]}), - Z(node=ancilla[1], domain={input_node}), - ) - ) - return ancilla[1], seq - - @classmethod - def _y_command(cls, input_node: int, ancilla: Sequence[int]) -> tuple[int, list[command.CommandType]]: - """MBQC commands for Pauli Y gate. - - Parameters - ---------- - input_node : int - input node index - ancilla : list of four ints - ancilla node indices to be added to graph - - Returns - ------- - out_node : int - control node on graph after the gate - commands : list - list of MBQC commands - """ - assert len(ancilla) == 4 - seq: list[CommandType] = [N(node=ancilla[0]), N(node=ancilla[1])] - seq.extend([N(node=ancilla[2]), N(node=ancilla[3])]) - seq.extend( - ( - E(nodes=(input_node, ancilla[0])), - E(nodes=(ancilla[0], ancilla[1])), - E(nodes=(ancilla[1], ancilla[2])), - E(nodes=(ancilla[2], ancilla[3])), - M(input_node, Measurement.Y), - M(ancilla[0], -Measurement.X, s_domain={input_node}), - M(ancilla[1], -Measurement.Y, s_domain={ancilla[0]}, t_domain={input_node}), - M(node=ancilla[2], s_domain={ancilla[1]}, t_domain={ancilla[0]}), - X(node=ancilla[3], domain={ancilla[2]}), - Z(node=ancilla[3], domain={ancilla[1]}), - ) - ) - return ancilla[3], seq - - @classmethod - def _z_command(cls, input_node: int, ancilla: Sequence[int]) -> tuple[int, list[command.CommandType]]: - """MBQC commands for Pauli Z gate. - - Parameters - ---------- - input_node : int - input node index - ancilla : list of two ints - ancilla node indices to be added to graph - - Returns - ------- - out_node : int - control node on graph after the gate - commands : list - list of MBQC commands - """ - assert len(ancilla) == 2 - seq: list[CommandType] = [N(node=ancilla[0]), N(node=ancilla[1])] - seq.extend( - ( - E(nodes=(input_node, ancilla[0])), - E(nodes=(ancilla[0], ancilla[1])), - M(input_node, -Measurement.X), - M(node=ancilla[0], s_domain={input_node}), - X(node=ancilla[1], domain={ancilla[0]}), - Z(node=ancilla[1], domain={input_node}), - ) - ) - return ancilla[1], seq - - @classmethod - def _rx_command( - cls, input_node: int, ancilla: Sequence[int], angle: ParameterizedAngle - ) -> tuple[int, list[command.CommandType]]: - """MBQC commands for X rotation gate. - - Parameters - ---------- - input_node : int - input node index - ancilla : list of two ints - ancilla node indices to be added to graph - angle : ParameterizedAngle - measurement angle in units of π - - Returns - ------- - out_node : int - control node on graph after the gate - commands : list - list of MBQC commands - """ - assert len(ancilla) == 2 - seq: list[CommandType] = [N(node=ancilla[0]), N(node=ancilla[1])] - seq.extend( - ( - E(nodes=(input_node, ancilla[0])), - E(nodes=(ancilla[0], ancilla[1])), - M(node=input_node), - M(ancilla[0], Measurement.XY(-angle), s_domain={input_node}), - X(node=ancilla[1], domain={ancilla[0]}), - Z(node=ancilla[1], domain={input_node}), - ) - ) - return ancilla[1], seq - - @classmethod - def _ry_command( - cls, input_node: int, ancilla: Sequence[int], angle: ParameterizedAngle - ) -> tuple[int, list[command.CommandType]]: - """MBQC commands for Y rotation gate. - - Parameters - ---------- - input_node : int - input node index - ancilla : list of four ints - ancilla node indices to be added to graph - angle : ParameterizedAngle - rotation angle in units of π - - Returns - ------- - out_node : int - control node on graph after the gate - commands : list - list of MBQC commands - """ - assert len(ancilla) == 4 - seq: list[CommandType] = [N(node=ancilla[0]), N(node=ancilla[1])] - seq.extend([N(node=ancilla[2]), N(node=ancilla[3])]) - seq.extend( - ( - E(nodes=(input_node, ancilla[0])), - E(nodes=(ancilla[0], ancilla[1])), - E(nodes=(ancilla[1], ancilla[2])), - E(nodes=(ancilla[2], ancilla[3])), - M(input_node, Measurement.Y), - M(ancilla[0], Measurement.XY(-angle), s_domain={input_node}), - M(ancilla[1], -Measurement.Y, s_domain={ancilla[0]}, t_domain={input_node}), - M(node=ancilla[2], s_domain={ancilla[1]}, t_domain={ancilla[0]}), - X(node=ancilla[3], domain={ancilla[2]}), - Z(node=ancilla[3], domain={ancilla[1]}), - ) - ) - return ancilla[3], seq - - @classmethod - def _rz_command( - cls, input_node: int, ancilla: Sequence[int], angle: ParameterizedAngle - ) -> tuple[int, list[command.CommandType]]: - """MBQC commands for Z rotation gate. - - Parameters - ---------- - input_node : int - input node index - ancilla : list of two ints - ancilla node indices to be added to graph - angle : ParameterizedAngle - measurement angle in units of π - - Returns - ------- - out_node : int - node on graph after the gate - commands : list - list of MBQC commands - """ - assert len(ancilla) == 2 - seq: list[CommandType] = [N(node=ancilla[0]), N(node=ancilla[1])] # assign new qubit labels - seq.extend( - ( - E(nodes=(input_node, ancilla[0])), - E(nodes=(ancilla[0], ancilla[1])), - M(input_node, Measurement.XY(-angle)), - M(node=ancilla[0], s_domain={input_node}), - X(node=ancilla[1], domain={ancilla[0]}), - Z(node=ancilla[1], domain={input_node}), - ) + indices: list[int | None] = list(range(self.width)) + n_nodes = self.width + measurements: dict[int, BlochMeasurement] = {} + classical_outputs: dict[int, command.M] = {} + inputs = list(range(n_nodes)) + graph: nx.Graph[int] = nx.Graph() + graph.add_nodes_from(inputs) + x_corrections: dict[int, set[int]] = {} + for instr in instructions_to_jcz(self.instruction): + match instr.kind: + case InstructionKind.M: + target = indices[instr.target] + if target is None: + raise RuntimeError("Ill-formed circuit") + classical_outputs[target] = command.M(target, PauliMeasurement(instr.axis)) + indices[instr.target] = None + continue + case InstructionKind.J: + target = indices[instr.target] + if target is None: + raise RuntimeError("Ill-formed circuit") + graph.add_edge(target, n_nodes) # Also adds nodes + measurements[target] = Measurement.XY(normalize_angle(-instr.angle)) + indices[instr.target] = n_nodes + x_corrections[target] = {n_nodes} # X correction on ancilla + n_nodes += 1 + continue + case InstructionKind.CZ: + t0, t1 = instr.targets + i0, i1 = indices[t0], indices[t1] + if i0 is None or i1 is None: + raise RuntimeError("Ill-formed circuit") + # If edge exists, remove it; else, add it + if graph.has_edge(i0, i1): + graph.remove_edge(i0, i1) + else: + graph.add_edge(i0, i1) + continue + case _: + assert_never(instr.kind) + outputs = [i for i in indices if i is not None] + outputs.extend(classical_outputs.keys()) # Necessary for flow-finding step + og = OpenGraph( + graph=graph, + input_nodes=inputs, + output_nodes=outputs, + measurements=measurements, ) - return ancilla[1], seq - - @classmethod - def _ccx_command( - cls, - control_node1: int, - control_node2: int, - target_node: int, - ancilla: Sequence[int], - ) -> tuple[int, int, int, list[command.CommandType]]: - """MBQC commands for CCX gate. + z_corrections: dict[int, set[int]] = {} + for node, correctors in x_corrections.items(): + z_targets = og.neighbors(correctors) - {node} + if z_targets: + z_corrections[node] = z_targets + partial_order_layers = _corrections_to_partial_order_layers(og, x_corrections, z_corrections) + f: CausalFlow[BlochMeasurement] = CausalFlow(og, x_corrections, partial_order_layers) + return TranspiledFlow(f, classical_outputs) + + def transpile(self, *, transpile_swaps: bool = True) -> TranspiledPattern: + """Transpile a circuit via J-∧z decomposition to a pattern. Parameters ---------- - control_node1 : int - first control node on graph - control_node2 : int - second control node on graph - target_node : int - target node on graph - ancilla : list of int - ancilla node indices to be added to graph + transpile_swaps: bool, optional + If ``True`` (the default), SWAP gates are eliminated by switching the qubits. + If ``False``, SWAP gates are transpiled into a sequence of three CNOT gates. Returns ------- - control_out1 : int - first control node on graph after the gate - control_out2 : int - second control node on graph after the gate - target_out : int - target node on graph after the gate - commands : list - list of MBQC commands + TranspiledPattern + The result of the transpilation: a pattern and classical outputs. """ - assert len(ancilla) == 18 - seq: list[CommandType] = [N(node=ancilla[i]) for i in range(18)] # assign new qubit labels - seq.extend( - ( - E(nodes=(target_node, ancilla[0])), - E(nodes=(ancilla[0], ancilla[1])), - E(nodes=(ancilla[1], ancilla[2])), - E(nodes=(ancilla[1], control_node2)), - E(nodes=(control_node1, ancilla[14])), - E(nodes=(ancilla[2], ancilla[3])), - E(nodes=(ancilla[14], ancilla[4])), - E(nodes=(ancilla[3], ancilla[5])), - E(nodes=(ancilla[3], ancilla[4])), - E(nodes=(ancilla[5], ancilla[6])), - E(nodes=(control_node2, ancilla[6])), - E(nodes=(control_node2, ancilla[9])), - E(nodes=(ancilla[6], ancilla[7])), - E(nodes=(ancilla[9], ancilla[4])), - E(nodes=(ancilla[9], ancilla[10])), - E(nodes=(ancilla[7], ancilla[8])), - E(nodes=(ancilla[10], ancilla[11])), - E(nodes=(ancilla[4], ancilla[8])), - E(nodes=(ancilla[4], ancilla[11])), - E(nodes=(ancilla[4], ancilla[16])), - E(nodes=(ancilla[8], ancilla[12])), - E(nodes=(ancilla[11], ancilla[15])), - E(nodes=(ancilla[12], ancilla[13])), - E(nodes=(ancilla[16], ancilla[17])), - M(node=target_node), - M(node=ancilla[0], s_domain={target_node}), - M(node=ancilla[1], s_domain={ancilla[0]}, t_domain={target_node}), - M(node=control_node1), - M(ancilla[2], Measurement.XY(-7 * ANGLE_PI / 4), s_domain={ancilla[1]}, t_domain={ancilla[0]}), - M(node=ancilla[14], s_domain={control_node1}), - M(node=ancilla[3], s_domain={ancilla[2]}, t_domain={ancilla[1], ancilla[14]}), - M(ancilla[5], Measurement.XY(-ANGLE_PI / 4), s_domain={ancilla[3]}, t_domain={ancilla[2]}), - M(control_node2, Measurement.XY(-ANGLE_PI / 4), t_domain={ancilla[5], ancilla[0]}), - M(node=ancilla[6], s_domain={ancilla[5]}, t_domain={ancilla[3]}), - M(node=ancilla[9], s_domain={control_node2}, t_domain={ancilla[14]}), - M(ancilla[7], Measurement.XY(-7 * ANGLE_PI / 4), s_domain={ancilla[6]}, t_domain={ancilla[5]}), - M(ancilla[10], Measurement.XY(-7 * ANGLE_PI / 4), s_domain={ancilla[9]}, t_domain={control_node2}), - M( - ancilla[4], - Measurement.XY(-ANGLE_PI / 4), - s_domain={ancilla[14]}, - t_domain={control_node1, control_node2, ancilla[2], ancilla[7], ancilla[10]}, - ), - M(node=ancilla[8], s_domain={ancilla[7]}, t_domain={ancilla[14], ancilla[6]}), - M(node=ancilla[11], s_domain={ancilla[10]}, t_domain={ancilla[9], ancilla[14]}), - M(ancilla[12], Measurement.XY(-ANGLE_PI / 4), s_domain={ancilla[8]}, t_domain={ancilla[7]}), - M( - node=ancilla[16], - s_domain={ancilla[4]}, - t_domain={ancilla[14]}, - ), - X(node=ancilla[17], domain={ancilla[16]}), - X(node=ancilla[15], domain={ancilla[11]}), - X(node=ancilla[13], domain={ancilla[12]}), - Z(node=ancilla[17], domain={ancilla[4]}), - Z(node=ancilla[15], domain={ancilla[10]}), - Z(node=ancilla[13], domain={ancilla[8]}), - ) - ) - return ancilla[17], ancilla[15], ancilla[13], seq + if not transpile_swaps: + return self.transpile_to_causal_flow().to_pattern() + swap = _transpile_swaps(self) + result = swap.circuit.transpile_to_causal_flow().to_pattern() + result.pattern.reorder_output_nodes(swap.swap_output_nodes(result.pattern.output_nodes)) + classical_outputs = swap.swap_classical_outputs(result.classical_outputs) + return TranspiledPattern(result.pattern, classical_outputs) @overload def simulate_statevector( @@ -1026,7 +595,7 @@ def simulate_statevector( else: _backend.add_nodes(range(self.width), input_state) - classical_measures = [] + classical_measures: list[Outcome] = [] for i in range(len(self.instruction)): instr = self.instruction[i] @@ -1064,6 +633,8 @@ def evolve(op: Matrix, qargs: Iterable[int]) -> None: evolve_single(Ops.ry(instr.angle), instr.target) case instruction.InstructionKind.RZ: evolve_single(Ops.rz(instr.angle), instr.target) + case instruction.InstructionKind.J: + evolve_single(Ops.j(instr.angle), instr.target) case instruction.InstructionKind.RZZ: evolve(Ops.rzz(instr.angle), [instr.control, instr.target]) case instruction.InstructionKind.CCX: @@ -1133,15 +704,255 @@ def transpile_measurements_to_z_axis(self) -> Circuit: circuit.add(instr) return circuit + def transpile_j_to_rzh(self) -> Circuit: + """Return an equivalent circuit where all J gates have been replaced with RZ and H gates.""" + new_circuit = Circuit(self.width) + for instr in self.instruction: + match instr.kind: + case InstructionKind.J: + new_circuit.add(instruction.RZ(target=instr.target, angle=instr.angle)) + new_circuit.add(instruction.H(target=instr.target)) + case _: + new_circuit.add(instr) + return new_circuit -def _transpile_rzz(instructions: Iterable[InstructionType]) -> Iterator[InstructionTypeWithoutRZZ]: - for instr in instructions: - if instr.kind == InstructionKind.RZZ: - yield instruction.CNOT(control=instr.control, target=instr.target) - yield instruction.RZ(target=instr.target, angle=instr.angle) - yield instruction.CNOT(control=instr.control, target=instr.target) - else: - yield instr + +def decompose_rzz(instr: instruction.RZZ) -> Iterator[instruction.CNOT | instruction.RZ]: + """Yield a decomposition of RZZ(α) gate as CNOT(control, target)·Rz(target, α)·CNOT(control, target). + + Parameters + ---------- + instr: the RZZ instruction to decompose. + + Returns + ------- + the decomposition. + + """ + yield instruction.CNOT(target=instr.target, control=instr.control) + yield instruction.RZ(instr.target, instr.angle) + yield instruction.CNOT(target=instr.target, control=instr.control) + + +def decompose_ccx( + instr: instruction.CCX, +) -> Iterator[instruction.H | instruction.CNOT | instruction.RZ]: + """Yield a decomposition of the CCX gate into H, CNOT, T and T-dagger gates. + + This decomposition of the Toffoli gate can be found in + Michael A. Nielsen and Isaac L. Chuang, + Quantum Computation and Quantum Information, + Cambridge University Press, 2000 + (p. 182 in the 10th Anniversary Edition). + + Parameters + ---------- + instr: the CCX instruction to decompose. + + Returns + ------- + the decomposition. + + """ + c0, c1, t = instr.controls[0], instr.controls[1], instr.target + yield instruction.H(t) + yield instruction.CNOT(control=c1, target=t) + yield instruction.RZ(t, -ANGLE_PI / 4) + yield instruction.CNOT(control=c0, target=t) + yield instruction.RZ(t, ANGLE_PI / 4) + yield instruction.CNOT(control=c1, target=t) + yield instruction.RZ(t, -ANGLE_PI / 4) + yield instruction.CNOT(control=c0, target=t) + yield instruction.RZ(c1, -ANGLE_PI / 4) + yield instruction.RZ(t, ANGLE_PI / 4) + yield instruction.CNOT(control=c0, target=c1) + yield instruction.H(t) + yield instruction.RZ(c1, -ANGLE_PI / 4) + yield instruction.CNOT(control=c0, target=c1) + yield instruction.RZ(c0, ANGLE_PI / 4) + yield instruction.RZ(c1, ANGLE_PI / 2) + + +def decompose_cnot(instr: instruction.CNOT) -> Iterator[instruction.H | instruction.CZ]: + """Yield a decomposition of the CNOT gate as H·∧z·H. + + Vincent Danos, Elham Kashefi, Prakash Panangaden, The Measurement Calculus, 2007. + + Parameters + ---------- + instr: the CNOT instruction to decompose. + + Returns + ------- + the decomposition. + + """ + yield instruction.H(instr.target) + yield instruction.CZ((instr.control, instr.target)) + yield instruction.H(instr.target) + + +def decompose_swap(instr: instruction.SWAP) -> Iterator[instruction.CNOT]: + """Yield a decomposition of the SWAP gate as CNOT(0, 1)·CNOT(1, 0)·CNOT(0, 1). + + Michael A. Nielsen and Isaac L. Chuang, + Quantum Computation and Quantum Information, + Cambridge University Press, 2000 + (p. 23 in the 10th Anniversary Edition). + + Parameters + ---------- + instr: the SWAP instruction to decompose. + + Returns + ------- + the decomposition. + + """ + yield instruction.CNOT(control=instr.targets[0], target=instr.targets[1]) + yield instruction.CNOT(control=instr.targets[1], target=instr.targets[0]) + yield instruction.CNOT(control=instr.targets[0], target=instr.targets[1]) + + +def decompose_y(instr: instruction.Y) -> Iterator[instruction.X | instruction.Z]: + """Return a decomposition of the Y gate as X·Z. + + Parameters + ---------- + instr: the Y instruction to decompose. + + Returns + ------- + the decomposition. + + """ + yield instruction.Z(instr.target) + yield instruction.X(instr.target) + + +def decompose_rx(instr: instruction.RX) -> Iterator[instruction.J]: + """Yield a J decomposition of the RX gate. + + The Rx(α) gate is decomposed into J(α)·H (that is to say, J(α)·J(0)). + Vincent Danos, Elham Kashefi, Prakash Panangaden, The Measurement Calculus, 2007. + + Parameters + ---------- + instr: the RX instruction to decompose. + + Returns + ------- + the decomposition. + + """ + yield instruction.J(instr.target, 0) + yield instruction.J(instr.target, instr.angle) + + +def decompose_ry(instr: instruction.RY) -> Iterator[instruction.J]: + """Yield a J decomposition of the RY gate. + + The Ry(α) gate is decomposed into J(0)·J(π/2)·J(α)·J(-π/2). + Vincent Danos, Elham Kashefi, Prakash Panangaden, Robust and parsimonious realisations of unitaries in the one-way + model, 2004. + + Parameters + ---------- + instr: the RY instruction to decompose. + + Returns + ------- + the decomposition. + + """ + yield instruction.J(target=instr.target, angle=-ANGLE_PI / 2) + yield instruction.J(target=instr.target, angle=instr.angle) + yield instruction.J(target=instr.target, angle=ANGLE_PI / 2) + yield instruction.J(target=instr.target, angle=0) + + +def decompose_rz(instr: instruction.RZ) -> Iterator[instruction.J]: + """Yield a J decomposition of the RZ gate. + + The Rz(α) gate is decomposed into H·J(α) (that is to say, J(0)·J(α)). + Vincent Danos, Elham Kashefi, Prakash Panangaden, The Measurement Calculus, 2007. + + Parameters + ---------- + instr: the RZ instruction to decompose. + + Returns + ------- + the decomposition. + + """ + yield instruction.J(target=instr.target, angle=instr.angle) + yield instruction.J(target=instr.target, angle=0) + + +def instructions_to_jcz(instrs: Iterable[InstructionType]) -> Iterator[instruction.J | instruction.CZ | instruction.M]: + """Yield a J-∧z decomposition of the instruction. + + Parameters + ---------- + instr: the instruction to decompose. + + Returns + ------- + the decomposition. + + """ + for instr in instrs: + match instr.kind: + case InstructionKind.J | InstructionKind.CZ | InstructionKind.M: + yield instr + case InstructionKind.I: + return + case InstructionKind.H: + yield instruction.J(instr.target, 0) + case InstructionKind.S: + yield from decompose_rz(instruction.RZ(instr.target, ANGLE_PI / 2)) + case InstructionKind.X: + yield from decompose_rx(instruction.RX(instr.target, ANGLE_PI)) + case InstructionKind.Y: + yield from instructions_to_jcz(decompose_y(instr)) + case InstructionKind.Z: + yield from decompose_rz(instruction.RZ(instr.target, ANGLE_PI)) + case InstructionKind.RX: + yield from decompose_rx(instr) + case InstructionKind.RY: + yield from decompose_ry(instr) + case InstructionKind.RZ: + yield from decompose_rz(instr) + case InstructionKind.CCX: + yield from instructions_to_jcz(decompose_ccx(instr)) + case InstructionKind.RZZ: + yield from instructions_to_jcz(decompose_rzz(instr)) + case InstructionKind.CNOT: + yield from instructions_to_jcz(decompose_cnot(instr)) + case InstructionKind.SWAP: + yield from instructions_to_jcz(decompose_swap(instr)) + case _: + assert_never(instr.kind) + + +def normalize_angle(angle: ParameterizedAngle) -> ParameterizedAngle: + r"""Return an equivalent angle in range :math:`[0, 2 \cdot \pi)` if ``angle`` is instantiated. + + Parameters + ---------- + angle: ParameterizedAngle + An angle. + + Returns + ------- + ParameterizedAngle + An equivalent angle in range :math:`[0, 2 \cdot \pi)` if ``angle`` is instantiated. + If ``angle`` is parameterized, ``angle`` is returned unchanged. + """ + if isinstance(angle, float): + return angle % (2 * ANGLE_PI) + return angle @dataclass(frozen=True) @@ -1151,24 +962,71 @@ class TranspileSwapsResult: circuit: Circuit """Circuit without SWAP gates.""" - qubits: tuple[int | None, ...] + outputs: tuple[OutputIndex, ...] """ Tuple which has the same width as the circuit and which for every qubit of the original circuit provides the index of the - corresponding qubit in the output of the returned - circuit. Measured qubits are mapped to ``None`` in this tuple. + corresponding qubit in the output of the swapped circuit + (either measured or not). """ + def extract_outputs(self, kind: OutputKind) -> tuple[int, ...]: + """Return the sequence of outputs of the given kind.""" + return tuple(output.index for output in self.outputs if output.kind == kind) + + def extract_output_node_indices(self) -> tuple[int, ...]: + """Return for each output node, sorted in the order of the original circuit, the index of the corresponding output node in the order of the swapped circuit.""" + reduced_index = {} + reduced_counter = 0 + for index, output in enumerate(self.outputs): + if output.kind == OutputKind.Qubit: + reduced_index[index] = reduced_counter + reduced_counter += 1 + return tuple(reduced_index[index] for index in self.extract_outputs(OutputKind.Qubit)) + + def swap_output_nodes(self, output_nodes: Sequence[Node]) -> tuple[Node, ...]: + """Reorder the output nodes of a pattern obtained from a swapped circuit to restore the qubit ordering of the original circuit.""" + return tuple(output_nodes[index] for index in self.extract_output_node_indices()) + + def swap_classical_outputs(self, classical_outputs: Sequence[Node]) -> tuple[int, ...]: + """Reorder the classical outpus of a pattern obtained from a swapped circuit to restore the output ordering of the original circuit.""" + return tuple(classical_outputs[index] for index in self.extract_outputs(OutputKind.Bit)) + + +class OutputKind(Enum): + """Specify whether a qubit is measured or not.""" + + Qubit = enum.auto() + Bit = enum.auto() + + +@dataclass(frozen=True) +class OutputIndex: + """Index of a swapped qubit. + + If the qubit is measured, ``kind`` equals to `OutputKind.Bit` and + ``index`` is the index of the measurement. + + If the qubit is not measured, ``kind`` equals to `OutputKind.qubit` + and ``index`` is the index of the qubit in the swapped circuit. + """ + + kind: OutputKind + index: int + class _TranspileSwapVisitor(InstructionVisitor): - qubits: list[int | None] + outputs: list[OutputIndex] def __init__(self, width: int) -> None: - self.qubits = list(range(width)) + self.outputs = [OutputIndex(OutputKind.Qubit, index) for index in range(width)] @override def visit_qubit(self, qubit: int) -> int: - return _check_target(self.qubits, qubit) + target = self.outputs[qubit] + if target.kind == OutputKind.Bit: + raise RuntimeError(f"Qubit {qubit} has already been measured.") + return target.index def transpile_swaps(circuit: Circuit) -> TranspileSwapsResult: @@ -1183,24 +1041,32 @@ def transpile_swaps(circuit: Circuit) -> TranspileSwapsResult: ------- TranspileSwapsResult The field ``circuit`` contains an equivalent circuit without - SWAP gates. - The field ``qubits`` contains a tuple which has the same width - as the circuit and which for every qubit of the original - circuit provides the index of the corresponding qubit in the - output of the returned circuit. Measured qubits are mapped to - ``None`` in this tuple. + SWAP gates. The field ``outputs`` contains a tuple which has + the same width as the circuit. For every qubit of the original + circuit, either the qubit is not measured, and ``outputs`` + provides the index of the corresponding qubit in the output of + the returned circuit; or the qubit has been measured, and + ``outputs`` provides the index of the measurement. """ new_circuit = Circuit(circuit.width) visitor = _TranspileSwapVisitor(circuit.width) + measurement_index = 0 for instr in circuit.instruction: if instr.kind == InstructionKind.SWAP: u, v = instr.targets - visitor.qubits[u], visitor.qubits[v] = visitor.qubits[v], visitor.qubits[u] + # We apply the visitor to check that the qubits have not been measured. + visitor.visit_qubit(u) + visitor.visit_qubit(v) + visitor.outputs[u], visitor.outputs[v] = visitor.outputs[v], visitor.outputs[u] else: new_circuit.add(instr.visit(visitor)) if instr.kind == InstructionKind.M: - visitor.qubits[instr.target] = None - return TranspileSwapsResult(new_circuit, tuple(visitor.qubits)) + visitor.outputs[instr.target] = OutputIndex(OutputKind.Bit, measurement_index) + measurement_index += 1 + return TranspileSwapsResult(new_circuit, tuple(visitor.outputs)) + + +_transpile_swaps = transpile_swaps @overload @@ -1219,20 +1085,20 @@ def _initialize_backend( @overload def _initialize_backend( - backend: DenseStateBackend[_DenseStateT_co], + backend: DenseStateBackend[_DenseStateT], branch_selector: BranchSelector | None, -) -> DenseStateBackend[_DenseStateT_co]: ... +) -> DenseStateBackend[_DenseStateT]: ... def _initialize_backend( - backend: DenseStateBackend[_DenseStateT_co] | _DenseStateBackendLiteral, + backend: DenseStateBackend[_DenseStateT] | _DenseStateBackendLiteral, branch_selector: BranchSelector | None, -) -> _BuiltinDenseStateBackend | DenseStateBackend[_DenseStateT_co]: +) -> _BuiltinDenseStateBackend | DenseStateBackend[_DenseStateT]: """Initialize backend for circuit simulation. Parameters ---------- - backend: :class:`graphix.sim.base_backend.DenseStateBackend[_DenseStateT_co]`, 'statevector', or 'densitymatrix' + backend: :class:`graphix.sim.base_backend.DenseStateBackend[_DenseStateT]`, 'statevector', or 'densitymatrix' Simulation backend branch_selector: :class:`BranchSelector` Branch selector used for measurements. Can only be specified if ``backend`` is not an already instantiated :class:`Backend` object. If ``None``, it defaults to :class:`RandomBranchSelector`. diff --git a/noxfile.py b/noxfile.py index 947c8e102..886811efa 100644 --- a/noxfile.py +++ b/noxfile.py @@ -113,19 +113,17 @@ class ReverseDependency: @nox.parametrize( "package", [ - ReverseDependency( - "https://github.com/thierry-martinez/graphix-stim-backend", branch="fix/graphix_498_remove_pauli" - ), + ReverseDependency("https://github.com/thierry-martinez/graphix-stim-backend", branch="fix/add_jcz"), ReverseDependency("https://github.com/TeamGraphix/graphix-symbolic"), ReverseDependency("https://github.com/TeamGraphix/graphix-qasm-parser"), ReverseDependency( "https://github.com/thierry-martinez/veriphix", doctest_modules=False, install_target=".[dev]", - branch="fix/graphix_498_remove_pauli", + branch="fix_reorder_and_refresh", ), ReverseDependency("https://github.com/TeamGraphix/graphix-ibmq", doctest_modules=False), - ReverseDependency("https://github.com/qat-inria/graphix-stim-compiler", branch="ps_dim"), + ReverseDependency("https://github.com/qat-inria/graphix-stim-compiler"), ], ) def tests_reverse_dependencies(session: Session, package: ReverseDependency) -> None: diff --git a/pyproject.toml b/pyproject.toml index 261a140a3..c5c3fcc41 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -195,8 +195,10 @@ extraPaths = ["./stubs"] [tool.coverage.report] exclude_also = [ "if TYPE_CHECKING:", - "raise NotImplementedError\\(.*\\)", + "case _:", + "raise NotImplementedError(\\(.*\\))?", "return NotImplemented", + "raise RuntimeError(\\(.*\\))", "typing_extensions.assert_never\\(.*\\)", "assert_never\\(.*\\)", "@abc.abstractmethod", diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index 7c67909c5..622c5e2a7 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -2,6 +2,7 @@ import copy import functools +import itertools import random from typing import TYPE_CHECKING @@ -15,13 +16,17 @@ from graphix.channels import KrausChannel, dephasing_channel, depolarising_channel from graphix.fundamentals import ANGLE_PI, Plane from graphix.ops import Ops +from graphix.random_objects import rand_state_vector from graphix.sim.density_matrix import DensityMatrix, DensityMatrixBackend from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, Statevec from graphix.simulator import DefaultMeasureMethod from graphix.states import BasicStates, PlanarState from graphix.transpiler import Circuit +from tests.test_statevec import permute_with_swap if TYPE_CHECKING: + from collections.abc import Sequence + from numpy.random import Generator from graphix.measurements import Outcome @@ -929,3 +934,21 @@ def test_measure(self, outcome: Outcome) -> None: else np.kron(np.array([[1, 0], [0, 0]]), np.ones((2, 2)) / 2) ) assert np.allclose(backend.state.rho, expected_matrix) + + +@pytest.mark.parametrize("permutation", itertools.permutations(range(3))) +def test_permute(fx_rng: Generator, permutation: Sequence[int]) -> None: + nqubits = len(permutation) + dm = DensityMatrix(rand_state_vector(nqubits, fx_rng)) + dm_ref = copy.copy(dm) + dm.permute(permutation) + permute_with_swap(dm_ref, permutation) + assert np.array_equal(dm.rho, dm_ref.rho) + + +def test_permute_bad_permutation() -> None: + dm = DensityMatrix(nqubit=2) + with pytest.raises(ValueError, match="Permutation has length"): + dm.permute([0]) + with pytest.raises(ValueError, match="not a permutation"): + dm.permute([1, 2]) diff --git a/tests/test_opengraph.py b/tests/test_opengraph.py index 64c2a45a8..2e62d4eb0 100644 --- a/tests/test_opengraph.py +++ b/tests/test_opengraph.py @@ -18,6 +18,7 @@ from graphix.fundamentals import ANGLE_PI, Axis, Plane from graphix.measurements import Measurement from graphix.opengraph import OpenGraph, OpenGraphError +from graphix.optimization import StandardizedPattern from graphix.parameter import Placeholder from graphix.pattern import Pattern from graphix.random_objects import rand_circuit @@ -1123,7 +1124,9 @@ def test_from_to_pattern(self, fx_rng: Generator) -> None: depth = 2 circuit = rand_circuit(n_qubits, depth, fx_rng) pattern_ref = circuit.transpile().pattern - pattern = pattern_ref.extract_opengraph().to_pattern() + pattern = StandardizedPattern.from_pattern( + pattern_ref.extract_opengraph().to_pattern() + ).to_space_optimal_pattern() for plane in {Plane.XY, Plane.XZ, Plane.YZ}: alpha = 2 * ANGLE_PI * fx_rng.random() diff --git a/tests/test_optimization.py b/tests/test_optimization.py index 42db33200..9d2e1305d 100644 --- a/tests/test_optimization.py +++ b/tests/test_optimization.py @@ -82,8 +82,10 @@ def test_remove_useless_domains(fx_bg: PCG64, jumps: int) -> None: pattern = circuit.transpile().pattern pattern.standardize() pattern.shift_signals() + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() pattern2 = remove_useless_domains(pattern) + pattern2 = StandardizedPattern.from_pattern(pattern2).to_space_optimal_pattern() state = pattern.simulate_pattern(rng=rng) state2 = pattern2.simulate_pattern(rng=rng) assert state.isclose(state2) diff --git a/tests/test_parameter.py b/tests/test_parameter.py index 88f68d07e..85d0e7ecf 100644 --- a/tests/test_parameter.py +++ b/tests/test_parameter.py @@ -170,6 +170,7 @@ def test_random_circuit_with_parameters(fx_bg: PCG64, jumps: int, use_xreplace: pattern = circuit.transpile().pattern pattern.standardize() pattern.shift_signals() + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() pattern.minimize_space() assignment: dict[Parameter, float] = {alpha: rng.uniform(high=2), beta: rng.uniform(high=2)} diff --git a/tests/test_pattern.py b/tests/test_pattern.py index 6a87c4016..f672f1ce6 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -19,6 +19,7 @@ from graphix.fundamentals import ANGLE_PI, Angle, Plane from graphix.measurements import BlochMeasurement, Measurement, Outcome, PauliMeasurement from graphix.opengraph import OpenGraph +from graphix.optimization import StandardizedPattern from graphix.pattern import Pattern, PatternError, RunnabilityError, RunnabilityErrorReason, shift_outcomes from graphix.random_objects import rand_circuit, rand_gate from graphix.sim.density_matrix import DensityMatrix @@ -70,10 +71,11 @@ def test_standardize(self, fx_rng: Generator) -> None: depth = 1 circuit = rand_circuit(nqubits, depth, fx_rng) pattern = circuit.transpile().pattern - + pattern = pattern.infer_pauli_measurements() pattern.standardize() assert pattern.is_standard() state = circuit.simulate_statevector().statevec + pattern.remove_pauli_measurements() state_mbqc = pattern.simulate_pattern(rng=fx_rng) assert state_mbqc.isclose(state) @@ -153,6 +155,7 @@ def test_minimize_space_with_gflow(self, fx_bg: PCG64, jumps: int) -> None: pattern = circuit.transpile().pattern pattern.standardize() pattern.shift_signals(method="mc") + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() pattern.minimize_space() state = circuit.simulate_statevector().statevec @@ -213,6 +216,7 @@ def test_shift_signals(self, fx_bg: PCG64, jumps: int) -> None: pattern.standardize() pattern.shift_signals(method="mc") assert pattern.is_standard() + pattern = StandardizedPattern.from_pattern(pattern).to_space_optimal_pattern() state = circuit.simulate_statevector().statevec state_mbqc = pattern.simulate_pattern(rng=rng) assert state_mbqc.isclose(state) @@ -231,6 +235,7 @@ def test_pauli_measurement_random_circuit( pattern = circuit.transpile().pattern pattern.standardize() pattern.shift_signals(method="mc") + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() pattern.minimize_space() state = circuit.simulate_statevector().statevec @@ -246,6 +251,7 @@ def test_pauli_measurement_random_circuit_all_paulis(self, fx_bg: PCG64, jumps: pattern = circuit.transpile().pattern pattern.standardize() pattern.shift_signals(method="mc") + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() input_node_set = set(pattern.input_nodes) assert not any( @@ -284,6 +290,7 @@ def test_pauli_measurement(self) -> None: pattern = circuit.transpile().pattern pattern.standardize() pattern.shift_signals(method="mc") + pattern = pattern.infer_pauli_measurements() pattern_opt = pattern.remove_pauli_measurements(copy=True) isolated_nodes = pattern_opt.extract_isolated_nodes() assert isolated_nodes == set() @@ -300,8 +307,9 @@ def test_pauli_measured_against_nonmeasured(self, fx_bg: PCG64, jumps: int) -> N depth = 2 circuit = rand_circuit(nqubits, depth, rng) pattern = circuit.transpile().pattern - pattern.standardize() + pattern.minimize_space() pattern1 = copy.deepcopy(pattern) + pattern1 = pattern1.infer_pauli_measurements() pattern1.remove_pauli_measurements() state = pattern.simulate_pattern(rng=rng) state1 = pattern1.simulate_pattern(rng=rng) @@ -407,6 +415,7 @@ def test_pauli_measurement_then_standardize(self, fx_bg: PCG64, jumps: int) -> N depth = 3 circuit = rand_circuit(nqubits, depth, rng) pattern = circuit.transpile().pattern + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() pattern.standardize() pattern.minimize_space() @@ -635,6 +644,7 @@ def test_compose_5(self, fx_rng: Generator) -> None: p2 = circuit_2.transpile().pattern # inputs: [0] p, _ = p1.compose(p2, mapping={0: 1, 1: 2, 2: 3}) + p = StandardizedPattern.from_pattern(p).to_space_optimal_pattern() circuit_12 = Circuit(1) circuit_12.h(0) @@ -673,6 +683,7 @@ def test_compose_7(self, fx_rng: Generator) -> None: circuit_1.rz(0, alpha) p1 = circuit_1.transpile().pattern p1.remove_input_nodes() + p1 = p1.infer_pauli_measurements() p1.remove_pauli_measurements() circuit_2 = Circuit(1) @@ -792,10 +803,12 @@ def test_extract_partial_order_layers_results(self) -> None: c = Circuit(1) c.rz(0, 0.2) p = c.transpile().pattern + p = p.infer_pauli_measurements() p.remove_pauli_measurements() assert p.extract_partial_order_layers() == (frozenset({1}), frozenset({0})) p = Pattern(cmds=[N(0), N(1), N(2), M(0), E((1, 2)), X(1, {0}), M(2, Measurement.XY(0.3))]) + p = p.infer_pauli_measurements() p.remove_pauli_measurements() assert p.extract_partial_order_layers() == (frozenset({1}), frozenset({2})) @@ -918,6 +931,8 @@ def test_extract_causal_flow_rnd_circuit(self, fx_bg: PCG64, jumps: int) -> None p_ref = circuit_1.transpile().pattern p_test = p_ref.to_bloch().extract_causal_flow().to_corrections().to_pattern().infer_pauli_measurements() + p_ref = p_ref.infer_pauli_measurements() + p_test = p_test.infer_pauli_measurements() p_ref.remove_pauli_measurements() p_test.remove_pauli_measurements() @@ -936,6 +951,8 @@ def test_extract_gflow_rnd_circuit(self, fx_bg: PCG64, jumps: int) -> None: p_ref = circuit_1.transpile().pattern p_test = p_ref.to_bloch().extract_gflow().to_corrections().to_pattern().infer_pauli_measurements() + p_ref = p_ref.infer_pauli_measurements() + p_test = p_test.infer_pauli_measurements() p_ref.remove_pauli_measurements() p_test.remove_pauli_measurements() @@ -1031,8 +1048,10 @@ def test_extract_xzc_rnd_circuit(self, fx_bg: PCG64, jumps: int) -> None: xzc.check_well_formed() p_test = xzc.to_pattern() - for p in [p_ref, p_test]: - p.remove_pauli_measurements() + p_ref = p_ref.infer_pauli_measurements() + p_ref.remove_pauli_measurements() + p_test = p_test.infer_pauli_measurements() + p_test.remove_pauli_measurements() s_ref = p_ref.simulate_pattern(rng=rng) s_test = p_test.simulate_pattern(rng=rng) diff --git a/tests/test_pretty_print.py b/tests/test_pretty_print.py index 28dedce55..6bf690d2b 100644 --- a/tests/test_pretty_print.py +++ b/tests/test_pretty_print.py @@ -119,7 +119,7 @@ def test_flow_pretty_print_random( flow_extractor: Callable[[OpenGraph[Measurement]], PauliFlow[Measurement]], ) -> None: rng = Generator(fx_bg.jumped(jumps)) - rand_og = rand_circuit(5, 5, rng=rng).transpile().pattern.extract_opengraph() + rand_og = rand_circuit(5, 5, rng=rng).transpile().pattern.infer_pauli_measurements().extract_opengraph() flow = flow_extractor(rand_og) flow.to_ascii() diff --git a/tests/test_qasm3_exporter.py b/tests/test_qasm3_exporter.py index e0b4d590d..651335e9c 100644 --- a/tests/test_qasm3_exporter.py +++ b/tests/test_qasm3_exporter.py @@ -55,7 +55,19 @@ def test_to_qasm3_random_circuit(fx_bg: PCG64, jumps: int) -> None: depth = 5 circuit = rand_circuit(nqubits, depth, rng=rng) pattern = circuit.transpile().pattern + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() pattern.minimize_space() - print(pattern) _qasm3 = pattern_to_qasm3(pattern) + + +def test_to_qasm3_failures() -> None: + circuit = Circuit(2) + circuit.m(0, Axis.X) + with pytest.raises(ValueError, match="OpenQASM3 only supports measurements on Z axis"): + circuit_to_qasm3(circuit, transpile=False) + circuit = circuit.transpile_measurements_to_z_axis() + circuit.j(1, 0.25) + with pytest.raises(ValueError, match="J gates must be decomposed before QASM3 export"): + circuit_to_qasm3(circuit, transpile=False) + _qasm3 = circuit_to_qasm3(circuit) diff --git a/tests/test_qasm3_exporter_to_graphix_parser.py b/tests/test_qasm3_exporter_to_graphix_parser.py index cd2421b06..9e31e9eb9 100644 --- a/tests/test_qasm3_exporter_to_graphix_parser.py +++ b/tests/test_qasm3_exporter_to_graphix_parser.py @@ -31,9 +31,10 @@ def check_round_trip(circuit: Circuit) -> None: qasm = circuit_to_qasm3(circuit) + check_circuit = circuit.transpile_j_to_rzh() parser = OpenQASMParser() parsed_circuit = parser.parse_str(qasm) - assert parsed_circuit.instruction == circuit.instruction + assert parsed_circuit.instruction == check_circuit.instruction @pytest.mark.parametrize("jumps", range(1, 11)) @@ -42,7 +43,7 @@ def test_circuit_to_qasm3(fx_bg: PCG64, jumps: int) -> None: nqubits = 5 depth = 4 # See https://github.com/TeamGraphix/graphix-qasm-parser/pull/5 - check_round_trip(rand_circuit(nqubits, depth, rng, use_cz=False)) + check_round_trip(rand_circuit(nqubits, depth, rng, use_j=True, use_cz=True)) @pytest.mark.parametrize( @@ -52,8 +53,7 @@ def test_circuit_to_qasm3(fx_bg: PCG64, jumps: int) -> None: instruction.RZZ(target=0, control=1, angle=ANGLE_PI / 4), instruction.CNOT(target=0, control=1), instruction.SWAP(targets=(0, 1)), - # See https://github.com/TeamGraphix/graphix-qasm-parser/pull/5 - # instruction.CZ(targets=(0, 1)), + instruction.CZ(targets=(0, 1)), instruction.H(target=0), instruction.S(target=0), instruction.X(target=0), @@ -63,7 +63,22 @@ def test_circuit_to_qasm3(fx_bg: PCG64, jumps: int) -> None: instruction.RX(target=0, angle=ANGLE_PI / 4), instruction.RY(target=0, angle=ANGLE_PI / 4), instruction.RZ(target=0, angle=ANGLE_PI / 4), + instruction.J(target=0, angle=ANGLE_PI / 4), ], ) def test_instruction_to_qasm3(instruction: InstructionType) -> None: check_round_trip(Circuit(3, instr=[instruction])) + + +def test_j_to_qasm3() -> None: + circuit = Circuit(3, instr=[instruction.J(target=0, angle=ANGLE_PI / 4)]) + qasm = circuit_to_qasm3(circuit) + parser = OpenQASMParser() + parsed_circuit = parser.parse_str(qasm) + assert parsed_circuit.instruction == circuit.transpile_j_to_rzh().instruction + + +def test_j_to_qasm3_failure() -> None: + circuit = Circuit(3, instr=[instruction.J(target=0, angle=ANGLE_PI / 4)]) + with pytest.raises(ValueError): + circuit_to_qasm3(circuit, transpile=False) diff --git a/tests/test_qasm3_exporter_to_qiskit.py b/tests/test_qasm3_exporter_to_qiskit.py index c44eb9dda..77c5bab76 100644 --- a/tests/test_qasm3_exporter_to_qiskit.py +++ b/tests/test_qasm3_exporter_to_qiskit.py @@ -117,8 +117,9 @@ def test_to_qasm3_random_circuit(fx_bg: PCG64, jumps: int) -> None: rng = Generator(fx_bg.jumped(jumps)) nqubits = 5 depth = 5 - circuit = rand_circuit(nqubits, depth, rng=rng) + circuit = rand_circuit(nqubits, depth, rng=rng, use_j=True) pattern = circuit.transpile().pattern + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() pattern.minimize_space() diff --git a/tests/test_remove_pauli_measurements.py b/tests/test_remove_pauli_measurements.py index b4729c488..10870ab36 100644 --- a/tests/test_remove_pauli_measurements.py +++ b/tests/test_remove_pauli_measurements.py @@ -257,6 +257,7 @@ def test_pattern_remove_pauli_measurements() -> None: circuit = Circuit(2) circuit.cnot(0, 1) pattern = circuit.transpile().pattern + pattern = pattern.infer_pauli_measurements() pattern2 = pattern.remove_pauli_measurements(copy=True) assert all_bloch_measurement_or_input_node( pattern2.input_nodes, (cmd for cmd in pattern2 if isinstance(cmd, Command.M)) diff --git a/tests/test_simulator.py b/tests/test_simulator.py index e10f03101..6c9ecbb81 100644 --- a/tests/test_simulator.py +++ b/tests/test_simulator.py @@ -2,36 +2,59 @@ from typing import TYPE_CHECKING -from graphix import Circuit -from graphix.sim.statevec import Statevec, StatevectorBackend -from graphix.states import BasicStates +import pytest + +from graphix import BasicStates, Pattern, Statevec, StatevectorBackend if TYPE_CHECKING: from numpy.random import Generator -def test_input_state_none(fx_rng: Generator) -> None: - circuit = Circuit(1) - circuit.h(0) - pattern = circuit.transpile().pattern - # By default, the initial state is |+>, therefore H|+> = |0>. - state = pattern.simulate_pattern(rng=fx_rng) +def test_no_explicit_input_state(hadamardpattern: Pattern, fx_rng: Generator) -> None: + # No explicit input state: the default initial state is |+⟩. + # H|+⟩ = |0⟩, so we expect the final state to be |0⟩. + state = hadamardpattern.simulate_pattern(rng=fx_rng) assert state.isclose(Statevec(BasicStates.ZERO)) - # With the initial state |0>, we obtain H|0> = |+>. - input_state = BasicStates.ZERO - state = pattern.simulate_pattern(input_state=input_state, rng=fx_rng) + + +def test_explicit_input_state_zero(hadamardpattern: Pattern, fx_rng: Generator) -> None: + # Provide an explicit input state |0⟩. + # H|0⟩ = |+⟩, so the final state should be |+⟩. + state = hadamardpattern.simulate_pattern(input_state=BasicStates.ZERO, rng=fx_rng) assert state.isclose(Statevec(BasicStates.PLUS)) - # With the initial state |0> prepared in the backend, we obtain - # H|0> = |+>. + + +def test_backend_prepared_zero(hadamardpattern: Pattern, fx_rng: Generator) -> None: + # Prepare the initial state in a backend and pass `input_state=None`. + # The backend already contains |0⟩ on its input nodes, + # therefore H|0⟩ = |+⟩. backend = StatevectorBackend() - backend.add_nodes(pattern.input_nodes, input_state) - state = pattern.simulate_pattern(backend=backend, input_state=None, rng=fx_rng) + backend.add_nodes(hadamardpattern.input_nodes, BasicStates.ZERO) + state = hadamardpattern.simulate_pattern(backend=backend, input_state=None, rng=fx_rng) assert state.isclose(Statevec(BasicStates.PLUS)) - # The backend already prepares |0>. If the simulator also prepares - # the input qubits in |+> (because we do not pass - # `input_state=None`), an additional qubit is introduced. The - # simulation therefore applies (I ⊗ H) on |0+>, resulting in |00>. + + +def test_no_prepared_qubits_and_input_state_none(hadamardpattern: Pattern, fx_rng: Generator) -> None: + # No prepared qubits in the backend and `input_state=None`. + # This is ambiguous, so a ValueError must be raised. + backend = StatevectorBackend() + with pytest.raises(ValueError, match="the backend is expected to have 1 input nodes already prepared"): + hadamardpattern.simulate_pattern(backend=backend, input_state=None, rng=fx_rng) + + +def test_prepared_qubits_and_input_state(hadamardpattern: Pattern, fx_rng: Generator) -> None: + # Backend already contains a state (|0⟩) **and** we ask the + # simulator to prepare its own input state (by omitting `input_state`). + # This would lead to double-allocation of qubits, so a ValueError is + # raised. + backend = StatevectorBackend() + backend.add_nodes(hadamardpattern.input_nodes, BasicStates.ZERO) + with pytest.raises(ValueError, match="the backend is expected to have no pre-allocated qubits"): + hadamardpattern.simulate_pattern(backend=backend, rng=fx_rng) + + +def test_node_index_after_finalize() -> None: + pattern = Pattern(input_nodes=[0, 1], output_nodes=[1, 0]) backend = StatevectorBackend() - backend.add_nodes(pattern.input_nodes, input_state) - state = pattern.simulate_pattern(backend=backend, rng=fx_rng) - assert state.isclose(Statevec(BasicStates.ZERO, nqubit=2)) + pattern.simulate_pattern(backend=backend) + assert list(backend.node_index) == [1, 0] diff --git a/tests/test_statevec.py b/tests/test_statevec.py index e7ce8d925..a91e56f7c 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -1,6 +1,8 @@ from __future__ import annotations +import copy import functools +import itertools from typing import TYPE_CHECKING import numpy as np @@ -8,15 +10,19 @@ from graphix.fundamentals import ANGLE_PI, Plane from graphix.pattern import Pattern +from graphix.random_objects import rand_state_vector +from graphix.sim.base_backend import NodeIndex from graphix.sim.statevec import Statevec, _norm_numeric from graphix.states import BasicStates, PlanarState if TYPE_CHECKING: - from collections.abc import Mapping + from collections.abc import Mapping, Sequence from typing import Literal from numpy.random import Generator + from graphix.sim.base_backend import DenseState + _ENCODING = Literal["LSB", "MSB"] @@ -243,3 +249,32 @@ def test_normalize() -> None: statevec = Statevec(nqubit=1, data=BasicStates.PLUS) statevec.remove_qubit(0) assert _norm_numeric(statevec.psi.astype(np.complex128, copy=False)) == 1 + + +def permute_with_swap(dense_state: DenseState, permutation: Sequence[int]) -> None: + nqubits = len(permutation) + node_index = NodeIndex() + node_index.extend(range(nqubits)) + for i, ind in enumerate(permutation): + if node_index.index(ind) != i: + move_from = node_index.index(ind) + dense_state.swap((i, move_from)) + node_index.swap(i, move_from) + + +@pytest.mark.parametrize("permutation", itertools.permutations(range(3))) +def test_permute(fx_rng: Generator, permutation: Sequence[int]) -> None: + nqubits = len(permutation) + statevec = Statevec(rand_state_vector(nqubits, fx_rng)) + statevec_ref = copy.copy(statevec) + statevec.permute(permutation) + permute_with_swap(statevec_ref, permutation) + assert np.array_equal(statevec.psi, statevec_ref.psi) + + +def test_permute_bad_permutation() -> None: + statevec = Statevec(nqubit=2) + with pytest.raises(ValueError, match="Permutation has length"): + statevec.permute([0]) + with pytest.raises(ValueError, match="not a permutation"): + statevec.permute([1, 2]) diff --git a/tests/test_tnsim.py b/tests/test_tnsim.py index c347e492b..2bb76a21b 100644 --- a/tests/test_tnsim.py +++ b/tests/test_tnsim.py @@ -330,6 +330,7 @@ def test_with_graphtrans(self, fx_bg: PCG64, jumps: int, fx_rng: Generator) -> N pattern = circuit.transpile().pattern pattern.standardize() pattern.shift_signals() + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() state = circuit.simulate_statevector().statevec tn_mbqc = pattern.simulate_pattern(backend="tensornetwork", rng=fx_rng) @@ -348,6 +349,7 @@ def test_with_graphtrans_sequential(self, fx_bg: PCG64, jumps: int, fx_rng: Gene pattern = circuit.transpile().pattern pattern.standardize() pattern.shift_signals() + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() state = circuit.simulate_statevector().statevec tn_mbqc = pattern.simulate_pattern(backend="tensornetwork", graph_prep="sequential", rng=fx_rng) @@ -396,6 +398,7 @@ def test_evolve(self, fx_bg: PCG64, jumps: int, fx_rng: Generator) -> None: pattern = circuit.transpile().pattern pattern.standardize() pattern.shift_signals() + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() state = circuit.simulate_statevector().statevec tn_mbqc = pattern.simulate_pattern(backend="tensornetwork", rng=fx_rng) diff --git a/tests/test_transpiler.py b/tests/test_transpiler.py index f848a5ed2..f3ae5f181 100644 --- a/tests/test_transpiler.py +++ b/tests/test_transpiler.py @@ -7,14 +7,15 @@ from numpy.random import PCG64, Generator from graphix import instruction -from graphix.branch_selector import ConstBranchSelector +from graphix.branch_selector import ConstBranchSelector, FixedBranchSelector from graphix.fundamentals import ANGLE_PI, Axis, Sign from graphix.instruction import I, InstructionKind from graphix.random_objects import rand_circuit, rand_gate, rand_state_vector from graphix.sim.density_matrix import DensityMatrix -from graphix.sim.statevec import Statevec +from graphix.sim.statevec import Statevec, StatevectorBackend +from graphix.simulator import DefaultMeasureMethod from graphix.states import BasicStates -from graphix.transpiler import Circuit, transpile_swaps +from graphix.transpiler import Circuit, OutputIndex, OutputKind, decompose_ccx, transpile_swaps from tests.test_branch_selector import CheckedBranchSelector if TYPE_CHECKING: @@ -43,111 +44,28 @@ lambda rng: instruction.RX(0, rng.random() * 2 * ANGLE_PI), lambda rng: instruction.RY(0, rng.random() * 2 * ANGLE_PI), lambda rng: instruction.RZ(0, rng.random() * 2 * ANGLE_PI), + lambda rng: instruction.J(0, rng.random() * 2 * ANGLE_PI), ] class TestTranspilerUnitGates: - def test_cz(self, fx_rng: Generator) -> None: - circuit = Circuit(2) - circuit.cz(0, 1) - pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) - - def test_cnot(self, fx_rng: Generator) -> None: - circuit = Circuit(2) - circuit.cnot(0, 1) - pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) - - def test_hadamard(self, fx_rng: Generator) -> None: - circuit = Circuit(1) - circuit.h(0) - pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) - - def test_s(self, fx_rng: Generator) -> None: - circuit = Circuit(1) - circuit.s(0) - pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) - - def test_x(self, fx_rng: Generator) -> None: - circuit = Circuit(1) - circuit.x(0) - pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) - - def test_y(self, fx_rng: Generator) -> None: - circuit = Circuit(1) - circuit.y(0) - pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) - - def test_z(self, fx_rng: Generator) -> None: - circuit = Circuit(1) - circuit.z(0) - pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) - - def test_rx(self, fx_rng: Generator) -> None: - theta = fx_rng.uniform() * 2 * ANGLE_PI - circuit = Circuit(1) - circuit.rx(0, theta) - pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) - - def test_ry(self, fx_rng: Generator) -> None: - theta = fx_rng.uniform() * 2 * ANGLE_PI - circuit = Circuit(1) - circuit.ry(0, theta) - pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) - - def test_rz(self, fx_rng: Generator) -> None: - theta = fx_rng.uniform() * 2 * ANGLE_PI - circuit = Circuit(1) - circuit.rz(0, theta) - pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) - - def test_i(self, fx_rng: Generator) -> None: - circuit = Circuit(1) - circuit.i(0) + @pytest.mark.parametrize("instruction", INSTRUCTION_TEST_CASES) + def test_instruction_flow(self, fx_rng: Generator, instruction: InstructionTestCase) -> None: + circuit = Circuit(3, instr=[instruction(fx_rng)]) pattern = circuit.transpile().pattern - state = circuit.simulate_statevector(rng=fx_rng).statevec - state_mbqc = pattern.simulate_pattern(rng=fx_rng) - assert state_mbqc.isclose(state) + circuit.transpile_to_causal_flow().flow.check_well_formed() + flow = pattern.to_bloch().extract_causal_flow() + flow.check_well_formed() @pytest.mark.parametrize("jumps", range(1, 11)) - def test_ccx(self, fx_bg: PCG64, jumps: int) -> None: + @pytest.mark.parametrize("instruction", INSTRUCTION_TEST_CASES) + def test_instructions(self, fx_bg: PCG64, jumps: int, instruction: InstructionTestCase) -> None: rng = Generator(fx_bg.jumped(jumps)) - nqubits = 4 - depth = 6 - circuit = rand_circuit(nqubits, depth, rng, use_ccx=True) + circuit = Circuit(3, instr=[instruction(rng)]) pattern = circuit.transpile().pattern - pattern.minimize_space() - state = circuit.simulate_statevector(rng=rng).statevec - state_mbqc = pattern.simulate_pattern(rng=rng) + input_state = rand_state_vector(3, rng=rng) + state = circuit.simulate_statevector(input_state=input_state).statevec + state_mbqc = pattern.simulate_pattern(input_state=input_state, rng=rng) assert state_mbqc.isclose(state) def test_transpiled(self, fx_rng: Generator) -> None: @@ -185,6 +103,21 @@ def test_measure( elif isinstance(state_mbqc, DensityMatrix) and isinstance(state, DensityMatrix): assert np.allclose(state_mbqc.rho, state.rho) + @pytest.mark.parametrize("jumps", range(1, 11)) + @pytest.mark.parametrize("axis", [Axis.X, Axis.Y, Axis.Z]) + @pytest.mark.parametrize("outcome", [0, 1]) + def test_measure_early(self, fx_bg: PCG64, jumps: int, axis: Axis, outcome: Outcome) -> None: + rng = Generator(fx_bg.jumped(jumps)) + circuit = Circuit(3) + circuit.m(0, axis) + circuit.cnot(1, 2) + input_state = rand_state_vector(3, rng=rng) + branch_selector = ConstBranchSelector(outcome) + state = circuit.simulate_statevector(rng=rng, input_state=input_state, branch_selector=branch_selector).statevec + pattern = circuit.transpile().pattern + state_mbqc = pattern.simulate_pattern(rng=rng, input_state=input_state, branch_selector=branch_selector) + assert state_mbqc.isclose(state) + @pytest.mark.parametrize("input_axis", [Axis.X, Axis.Y, Axis.Z]) @pytest.mark.parametrize("input_sign", [Sign.PLUS, Sign.MINUS]) @pytest.mark.parametrize("measurement_axis", [Axis.X, Axis.Y, Axis.Z]) @@ -227,6 +160,130 @@ def test_transpile_measurements_to_z_axis(self, fx_bg: PCG64, jumps: int, axis: ).statevec assert state_z.isclose(state) + @pytest.mark.parametrize("jumps", range(1, 11)) + def test_transpile_j_to_rzh(self, fx_bg: PCG64, jumps: int) -> None: + rng = Generator(fx_bg.jumped(jumps)) + nqubits = 3 + depth = 2 + circuit = rand_circuit(nqubits, depth, rng, use_j=True, use_ccx=True, use_rzz=True) + circuit.j(0, 0.5) # Ensure that there is at least one J instruction + assert any(instr.kind == InstructionKind.J for instr in circuit.instruction) + circuit2 = circuit.transpile_j_to_rzh() + assert not any(instr.kind == InstructionKind.J for instr in circuit2.instruction) + state = circuit.simulate_statevector(rng=rng).statevec + state2 = circuit2.simulate_statevector(rng=rng).statevec + assert state.fidelity(state2) == pytest.approx(1) + + @pytest.mark.parametrize("jumps", range(1, 11)) + @pytest.mark.parametrize("axis", [Axis.X, Axis.Y, Axis.Z]) + @pytest.mark.parametrize("outcome", [0, 1]) + def test_transpile_swaps_with_measurements(self, fx_bg: PCG64, jumps: int, axis: Axis, outcome: Outcome) -> None: + rng = Generator(fx_bg.jumped(jumps)) + circuit = Circuit(3) + circuit.swap(0, 1) + circuit.swap(0, 2) + circuit.cnot(1, 2) + circuit.m(1, axis) + circuit.i(0) + transpiled_swaps = transpile_swaps(circuit) + circuit2 = transpiled_swaps.circuit + assert not any(instr.kind == InstructionKind.SWAP for instr in circuit2.instruction) + assert I(2) in circuit2.instruction + input_state = rand_state_vector(3, rng=rng) + branch_selector = ConstBranchSelector(outcome) + state = circuit.simulate_statevector(rng=rng, input_state=input_state, branch_selector=branch_selector).statevec + state2 = circuit2.simulate_statevector( + rng=rng, input_state=input_state, branch_selector=branch_selector + ).statevec + assert transpiled_swaps.outputs == ( + OutputIndex(OutputKind.Qubit, 2), + OutputIndex(OutputKind.Bit, 0), + OutputIndex(OutputKind.Qubit, 1), + ) + state2.swap((0, 1)) + assert state.isclose(state2) + + def test_cz_ccx(self, fx_rng: Generator) -> None: + """Test case reported in issue #2. + + https://github.com/qat-inria/graphix-jcz-transpiler/issues/2 + """ + circuit = Circuit(width=3) + circuit.cz(2, 0) + circuit.ccx(0, 1, 2) + ref_state = circuit.simulate_statevector(rng=fx_rng).statevec + pattern = circuit.transpile().pattern + state = pattern.simulate_pattern(rng=fx_rng) + assert state.isclose(ref_state) + + def test_ccx_decomposition(self) -> None: + circuit = Circuit(width=3) + circuit.cz(2, 0) + circuit.ccx(0, 1, 2) + circuit2 = Circuit(width=3) + circuit2.cz(2, 0) + circuit2.extend(decompose_ccx(instruction.CCX(controls=(0, 1), target=2))) + state = circuit.simulate_statevector().statevec + state2 = circuit2.simulate_statevector().statevec + assert state.isclose(state2) + + def test_cnot_cz(self, fx_rng: Generator) -> None: + """Test regression about output node reordering.""" + circuit = Circuit(width=3, instr=[instruction.CNOT(0, 1), instruction.CZ((0, 1))]) + state = circuit.simulate_statevector(rng=fx_rng).statevec + pattern = circuit.transpile().pattern + state_mbqc = pattern.simulate_pattern(rng=fx_rng) + assert state.isclose(state_mbqc) + + @pytest.mark.parametrize("jumps", range(1, 6)) + @pytest.mark.parametrize("axes", [[Axis.X, Axis.Y], [Axis.X, Axis.Y, Axis.Z]]) + def test_classical_outputs_consistency(self, fx_bg: PCG64, jumps: int, axes: list[Axis]) -> None: + """Check that `classical_outputs` are in the same order as `classical_measures`.""" + rng = Generator(fx_bg.jumped(jumps)) + n = len(axes) + width = n + 1 + circuit = Circuit(width) + for q in range(n): + circuit.cnot(q, q + 1) + for q, axis in enumerate(axes): + circuit.m(q, axis) + + transpile_result = circuit.transpile() + pattern = transpile_result.pattern + expected_outcomes: list[Outcome] = [1 if q % 2 else 0 for q in range(n)] + results_circuit: dict[int, Outcome] = dict(zip(range(n), expected_outcomes, strict=False)) + m_outcomes = dict(zip(transpile_result.classical_outputs, expected_outcomes, strict=False)) + non_output_nodes = pattern.extract_nodes() - set(pattern.output_nodes) + results_pattern: dict[int, Outcome] = {node: m_outcomes.get(node, 0) for node in non_output_nodes} + input_state = rand_state_vector(width, rng=rng) + measure_method = DefaultMeasureMethod() + circuit_result = circuit.simulate_statevector( + rng=rng, + input_state=input_state, + branch_selector=FixedBranchSelector(results=results_circuit), + ) + pattern.simulate_pattern( + rng=rng, + input_state=input_state, + branch_selector=FixedBranchSelector(results=results_pattern), + measure_method=measure_method, + ) + assert len(transpile_result.classical_outputs) == len(circuit_result.classical_measures) + pattern_measures = [measure_method.results[node] for node in transpile_result.classical_outputs] + assert pattern_measures == list(circuit_result.classical_measures) + assert pattern_measures == expected_outcomes + + def test_classical_outputs_empty(self) -> None: + """Circuits with no M instructions produce empty classical_outputs.""" + circuit = Circuit(2) + circuit.cnot(0, 1) + circuit.h(0) + result = circuit.transpile() + assert len(result.classical_outputs) == 0 + assert len(circuit.simulate_statevector().classical_measures) == 0 + + +class TestCircuits: def test_add_extend(self) -> None: circuit = Circuit(3) circuit.ccx(0, 1, 2) @@ -300,11 +357,7 @@ def test_transpile_swaps(fx_bg: PCG64, jumps: int) -> None: assert not any(instr.kind == InstructionKind.SWAP for instr in circuit2.instruction) state = circuit.simulate_statevector(rng=rng).statevec state2 = circuit2.simulate_statevector(rng=rng).statevec - qubits: list[int] = [] - for qubit in transpiled_swaps.qubits: - assert qubit is not None - qubits.append(qubit) - state2.psi = np.transpose(state2.psi, qubits) + state2.permute(transpiled_swaps.extract_output_node_indices()) assert state.isclose(state2) @@ -327,6 +380,37 @@ def test_transpile_swaps_with_measurements(fx_bg: PCG64, jumps: int, axis: Axis, branch_selector = ConstBranchSelector(outcome) state = circuit.simulate_statevector(rng=rng, input_state=input_state, branch_selector=branch_selector).statevec state2 = circuit2.simulate_statevector(rng=rng, input_state=input_state, branch_selector=branch_selector).statevec - assert transpiled_swaps.qubits == (2, None, 1) + assert transpiled_swaps.outputs == ( + OutputIndex(OutputKind.Qubit, 2), + OutputIndex(OutputKind.Bit, 0), + OutputIndex(OutputKind.Qubit, 1), + ) state2.swap((0, 1)) assert state.isclose(state2) + + +def test_transpile_double_cz() -> None: + circuit = Circuit(2) + circuit.cz(0, 1) + circuit.cz(1, 0) + cf = circuit.transpile_to_causal_flow() + assert len(cf.flow.og.graph.edges) == 0 + + +def test_transpile_swaps_vs_no_transpile_swaps() -> None: + circuit = Circuit(2) + circuit.rx(0, 0.25) + circuit.ry(0, 0.25) + circuit.cz(0, 1) + circuit.swap(0, 1) + pattern_without_swap = circuit.transpile().pattern + pattern_with_swap = circuit.transpile(transpile_swaps=False).pattern + state_without_swap = pattern_without_swap.simulate_pattern() + state_with_swap = pattern_with_swap.simulate_pattern() + assert state_without_swap.isclose(state_with_swap) + + +def test_backend_branch_selector() -> None: + circ = Circuit(1) + with pytest.raises(ValueError, match="already instantiated"): + circ.simulate_statevector(backend=StatevectorBackend(), branch_selector=ConstBranchSelector(0)) diff --git a/tests/test_visualization.py b/tests/test_visualization.py index c2f6942b2..598082bc8 100644 --- a/tests/test_visualization.py +++ b/tests/test_visualization.py @@ -158,6 +158,7 @@ def example_hadamard() -> Pattern: def example_local_clifford() -> Pattern: pattern = example_hadamard() + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() return pattern @@ -257,6 +258,7 @@ def test_draw_graph_reference(flow_and_not_pauli_presimulate: bool) -> Figure: # to have causal flow. pattern = pattern.to_bloch() else: + pattern = pattern.infer_pauli_measurements() pattern.remove_pauli_measurements() pattern.standardize() pattern.draw(