Encoding circuits

docs/Encoders.md

@@ -13,16 +13,15 @@ mystnb:
 
 # Encoder Circuit Synthesis for CSS Codes
 
-QECC provides functionality for synthesizing encoding circuits of arbitrary CSS codes. An encoder for an $[[n,k,d]]$ code is an isometry that encodes $k$ logical qubits into $n$ physical qubits.
+QECC provides functionality for synthesizing encoding circuits of arbitrary Stabilizer codes. An encoder for an $[[n,k,d]]$ code is an isometry that encodes $k$ logical qubits into $n$ physical qubits.
 
 Let's consider the synthesis of the encoding circuit of the $[[7,1,3]]$ Steane code.
 
 ```{code-cell} ipython3
 from mqt.qecc import CSSCode
 from mqt.qecc.circuit_synthesis import (
     depth_optimal_encoding_circuit,
-    gate_optimal_encoding_circuit,
-    heuristic_encoding_circuit,
+    gate_optimal_encoding_circuit
 )
 
 steane_code = CSSCode.from_code_name("steane")
@@ -31,14 +30,15 @@ print("Stabilizers:\n")
 print(steane_code.stabs_as_pauli_strings())
 print("\nLogicals:\n")
 print(steane_code.x_logicals_as_pauli_strings())
+print(steane_code.z_logicals_as_pauli_strings())
 ```
 
 There is not a unique encoding circuit but usually we would like to obtain an encoding circuit that is optimal with respect to some metric. QECC has functionality for synthesizing gate- or depth-optimal encoding circuits.
 
 Under the hood, this uses the SMT solver [z3](https://github.com/Z3Prover/z3). Of course this method scales only up to a few qubits. Synthesizing depth-optimal circuits is usually faster than synthesizing gate-optimal circuits.
 
 ```{code-cell} ipython3
-depth_opt = depth_optimal_encoding_circuit(steane_code, max_timeout=5)
+depth_opt = depth_optimal_encoding_circuit(steane_code, max_timeout=2)
 q_enc = depth_opt.get_uninitialized()
 
 print(f"Encoding qubits are qubits {q_enc}.")
@@ -49,7 +49,7 @@ depth_opt.draw('mpl')
 ```
 
 ```{code-cell} ipython3
-gate_opt = gate_optimal_encoding_circuit(steane_code, max_timeout=5)
+gate_opt = gate_optimal_encoding_circuit(steane_code, max_timeout=2)
 q_enc = gate_opt.get_uninitialized()
 
 print(f"Encoding qubits are qubits {q_enc}.")
@@ -66,106 +66,129 @@ In addition to the circuit, the synthesis methods also return the encoding qubit
 For larger codes, synthesizing optimal circuits is not feasible. In this case, QECC provides a heuristic synthesis method that tries to use as few CNOTs with the lowest depth as possible.
 
 ```{code-cell} ipython3
-heuristic_circ = heuristic_encoding_circuit(steane_code)
-q_enc = heuristic_circ.get_uninitialized()
+from mqt.qecc.circuit_synthesis import (
+    synthesize_encoding_circuit,
+)
 
-print(f"Encoding qubits are qubits {q_enc}.")
+heuristic_circ = synthesize_encoding_circuit(steane_code)
+q_enc = heuristic_circ.inputs()
+
+print(f"Messaging (logical input) qubits: {q_enc}")
 print(f"Circuit has depth {heuristic_circ.depth()}.")
 print(f"Circuit has {heuristic_circ.num_cnots()} CNOTs.")
 
 heuristic_circ.draw('mpl')
 ```
 
-## Synthesizing Encoders for Concatenated Codes
+By default the heuristic synthesis tries to optimize for two-qubit gate count. We can also tell the synthesis to optimize for depth.
 
-Encoders for concatenated codes can be constructed by concatenating encoding circuits. We can concatenate the $[[4,2,2]]$ code (with stabilizer generators $XXXX$ and $ZZZZ$) with itself by encoding $4$ qubits into two blocks of the code and then encoding these qubits one more time. This gives an $[[8,4,2]]$ code. The distance is still $2$ but if done the right way, some minimal-weight logicals have weight $4$.
+```{code-cell} ipython3
+from mqt.qecc.circuit_synthesis import (
+    SynthesisConfig
+)
 
-As an exercise, let's construct the concatenated circuit.
+config = SynthesisConfig(optimization_criterion="depth")
+heuristic_circ = synthesize_encoding_circuit(steane_code)
+q_enc = heuristic_circ.inputs()
 
-We start off by defining the code:
+print(f"Messaging (logical input) qubits: {q_enc}")
+print(f"Circuit has depth {heuristic_circ.depth()}.")
+print(f"Circuit has {heuristic_circ.num_cnots()} CNOTs.")
 
-```{code-cell} ipython3
-import numpy as np
+heuristic_circ.draw('mpl')
+```
 
-d = 2
-x_stabs = np.ones((1, 4), dtype=np.int8)
-z_stabs = x_stabs
-code = CSSCode(x_stabs, z_stabs, d)
+The `inputs()` method returns a list of physical qubit indices representing the encoded logical information. All other qubits are ancillas initialized in the $|0\rangle$ or $|+\rangle$ state.
 
-print("Stabilizers:\n")
-print(code.stabs_as_pauli_strings())
-print("\nLogicals:\n")
-print(code.x_logicals_as_pauli_strings())
-print(code.z_logicals_as_pauli_strings())
-```
+# Encoder Circuit Synthesis for non-CSS Stabilizer Codes
 
-We have to be careful with the logicals. Each _anticommuting_ pair of logicals defines one logical qubit.
+QECC also supports encoding circuit synthesis for non-CSS stabilizer codes. For these codes, encoding circuits are built from two-qubit transvections (see arXiv:2503.14660 for details).
 
-As before, we synthesize the encoding circuit:
+Let's consider the $[[5,1,3]]$ code.
 
 ```{code-cell} ipython3
-encoder = depth_optimal_encoding_circuit(code, max_timeout=5)
-q_enc = encoder.get_uninitialized()
+from mqt.qecc import StabilizerCode
+from mqt.qecc.circuit_synthesis import synthesize_encoding_circuit
+
+stabs = ["XZZXI", "IXZZX", "XIXZZ", "ZXIXZ"]
+x_logicals = ["XXXXX"]
+z_logicals = ["ZZZZZ"]
+five_qubit_code = StabilizerCode(stabs, x_logicals=x_logicals, z_logicals=z_logicals)
+
+print("Stabilizers:")
+for stab in stabs:
+    print(f"  {stab}")
+print("\nX Logical:")
+print(f"  {x_logicals[0]}")
+print("Z Logical:")
+print(f"  {z_logicals[0]}")
+```
 
-print(f"Encoding qubits are qubits {q_enc}.")
-print(f"Circuit has depth {encoder.depth()}.")
-print(f"Circuit has {encoder.num_cnots()} CNOTs.")
+The same `synthesize_encoding_circuit` function works for non-CSS codes. This method returns a `CliffordIsometry` object.
 
-encoder.draw('mpl')
+```{code-cell} ipython3
+encoder = synthesize_encoding_circuit(five_qubit_code)
+message_qubits = encoder.inputs()
+
+print(f"Message qubits (logical to physical mapping): {message_qubits}")
+print(f"Circuit has two-qubit depth {encoder.depth()}.")
+print(f"Circuit has {encoder.num_two_qubit_gates()} two-qubit gates.")
+
+encoder.draw(output='mpl', fold=False)
 ```
 
-Propagating Paulis from the encoding qubits at the input to the output will not necessarily yield the exact logicals given above. But the logical operators will be stabilizer equivalent.
+For displaying the circuit, the transvections are decomposed into CZ and single-qubit Clifford gates.
 
-Concatenating the circuits can be done as follows with qiskit:
+For non-CSS codes, depth-optimal synthesis is also available:
 
 ```{code-cell} ipython3
-from mqt.qecc.circuit_synthesis.circuits import compose_cnot_circuits
+from mqt.qecc.circuit_synthesis import depth_optimal_encoding_circuit_non_css
 
-n = 4
+for d in range(3, 9):
+    result = depth_optimal_encoding_circuit_non_css(five_qubit_code, max_depth=d, exact_two_qubit_count=True, max_two_qubit_gates=6)
 
-first_layer = encoder
-second_layer, mapping1, mapping2 = compose_cnot_circuits(encoder, encoder) # vertically composes circuits
+    if result == "UNSAT":
+        print(f"No solution for max_depth={d}")
+    else:
+        encoder = result
+        break
+print(f"Message qubits (logical to physical mapping): {message_qubits}")
+print(f"Circuit has two-qubit depth {encoder.depth()}.")
+print(f"Circuit has {encoder.num_two_qubit_gates()} two-qubit gates.")
 
-wiring = {0: mapping1[q_enc[0]], 1: mapping1[q_enc[1]], 2: mapping2[q_enc[0]], 3: mapping2[q_enc[1]]}
-concatenated, _, _ = compose_cnot_circuits(first_layer, second_layer, wiring)
-
-q_enc = concatenated.get_uninitialized()
-print(f"Encoding qubits are qubits {q_enc}.")
-print(f"Circuit has depth {concatenated.depth()}.")
-print(f"Circuit has {concatenated.num_cnots()} CNOTs.")
-
-concatenated.draw('mpl')
+encoder.draw(output='mpl', fold=False)
 ```
 
-Qubits $1$ and $2$ are still the encoding qubits and if we propagate Pauli $X$ and $Z$ to the output, we find that this is indeed the encoder for an $[[8,2,2]]$ code.
+This method uses SMT-based synthesis to find a depth-optimal encoding circuit, similar to the CSS case. The `max_depth` parameter limits the search depth. If no solution is found, it returns `"UNSAT"`.
+
+## Tweaking Parameters for Heuristic Synthesis
 
-This circuit has $3$ times as many CNOT gates as the encoder for the unconcatenated code because we needed to encode 3 times. Instead of concatenating the encoder circuits we can synthesize the encoders directly from the stabilizers of the concatenated code.
-We can obtain the code defined by the circuit directly from the circuit object.
+Let's consider a slightly larger example, the $[[23,1,7]]$ [Golay code](https://errorcorrectionzoo.org/c/qubit_golay).
 
 ```{code-cell} ipython3
-concatenated_code = concatenated.get_code()
+code = CSSCode.from_code_name("golay")
 
-print("Stabilizers:\n")
-print(concatenated_code.stabs_as_pauli_strings())
+encoder = synthesize_encoding_circuit(code)
+print(f"Messaging (logical input) qubits: {encoder.inputs()}")
+print(f"Circuit has depth {encoder.depth()}.")
+print(f"Circuit has {encoder.num_cnots()} CNOTs.")
 
-print("\nLogicals:\n")
-print(concatenated_code.x_logicals_as_pauli_strings())
-print(concatenated_code.z_logicals_as_pauli_strings())
+encoder.draw(output='mpl', fold=False, scale=0.5)
 ```
 
-Now we can directly synthesize the encoder:
+The way the greedy synthesis works in QECC is by trying to reduce the check matrix (or stabilizer tableau) of the code in question using as few elementary matrix operations (gates) as possible. The synthesis is greedily guided by some metric computed on the check matrix - the number of remaining entries in the check matrix, for example. Since the synthesis algorithm makes local choices, the search might go down branches of the search-tree leading to sub-optimal solution. QECC also provides a costlier search procedure which tries to look ahead which candidates in the search will lead to good results by completing the entire greedy synthesis for these candidates and making a choice based on the resulting circuits (as opposed to using the check matrix as a proxy for estimating).
+
+This search is generally costlier, but can lead to significantly better results. We can tell the synthesis to perform rollout-based synthesis by setting the appropriate flags in the config. `rollout` is an int parameter that determines for how many layers the search should perform the rollout. If it is set to `0`, no rollout is performed (default). If it is set to `1`, the `num_rollout_candidates` parameter determines for how many candidates per gate in the search the rollout is performed. This determines how many complete circuits are synthesized before the single gate is chosen leading to the locally best circuit. If `enable_early_termination` rollout is only performed until no better solutions are found. In that case, the search returns whatever the current best circuit is. If it is set to `False`, rollout will be performed until the last gate of the search is placed. This will take longer but leads to better results in general:
 
 ```{code-cell} ipython3
-encoder_concat_direct = depth_optimal_encoding_circuit(
-    concatenated_code, max_timeout=5
-)
-q_enc = encoder_concat_direct.get_uninitialized()
+from mqt.qecc.circuit_synthesis import SynthesisConfig
 
-print(f"Encoding qubits are qubits {q_enc}.")
-print(f"Circuit has depth {encoder_concat_direct.depth()}.")
-print(f"Circuit has {encoder_concat_direct.num_cnots()} CNOTs.")
+config = SynthesisConfig(rollout=1, num_rollout_candidates=5, optimization_criterion="gates", enable_early_termination=False)
 
-encoder_concat_direct.draw('mpl')
-```
+encoder = synthesize_encoding_circuit(code, config=config)
+print(f"Messaging (logical input) qubits: {encoder.inputs()}")
+print(f"Circuit has depth {encoder.depth()}.")
+print(f"Circuit has {encoder.num_cnots()} CNOTs.")
 
-We see that the circuit is more compact then the naively concatenated one. This is because the synthesis method exploits redundancy in the check matrix of the concatenated code.
+encoder.draw(output='mpl', fold=False, scale=0.5)
+```

notebooks/Correction Network.ipynb

@@ -0,0 +1,86 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "9f21998a-6ce6-4d08-a967-6cb687b9dea0",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from mqt.qecc.circuit_synthesis import correction_network"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "1ce0e5d2-d7dc-414a-ad41-fae660d866fc",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# remap qubits to be contiguous and store inverse maps\n",
+    "flip_index_map = {2: 0, 7: 1, 11: 2, 14: 3}\n",
+    "flip_index_map_inverse = {i: j for j, i in flip_index_map.items()}\n",
+    "correction_index_map = {0: 0, 3: 1, 4: 2, 6: 3, 8: 4, 10: 5, 12: 6}\n",
+    "correction_index_map_inverse = {i: j for j, i in correction_index_map.items()}\n",
+    "\n",
+    "# map from flipped qubits to wanted corrections\n",
+    "correction_map = {\n",
+    "    (11,): [0, 3, 4, 6, 8, 10, 12],\n",
+    "    (2,): [3, 4, 6, 10],\n",
+    "    (14,): [0, 3, 12],\n",
+    "    (7,): [3, 8],\n",
+    "    (11, 14): [4],\n",
+    "    (2, 14): [0, 10],\n",
+    "    (7, 11): [0, 6],\n",
+    "    (2, 7): [0, 4, 8],\n",
+    "    (7, 14): [4, 12],\n",
+    "    (2, 11, 14): [6, 12],\n",
+    "    (2, 7, 11): [10],\n",
+    "    (7, 11, 14): [6, 8],\n",
+    "    (2, 7, 14): [8, 10, 12],\n",
+    "    (2, 7, 11, 14): [3],\n",
+    "    (2, 11): [],\n",
+    "}\n",
+    "\n",
+    "correction_map_mapped = {}\n",
+    "for pattern, correction in correction_map.items():\n",
+    "    pattern_mapped = tuple(flip_index_map[q] for q in pattern)\n",
+    "    correction_mapped = [correction_index_map[q] for q in correction]\n",
+    "    correction_map_mapped[pattern_mapped] = correction_mapped\n",
+    "\n",
+    "n_data = 7\n",
+    "n_flips = 4\n",
+    "\n",
+    "network = correction_network(n_data, n_flips, correction_map_mapped)\n",
+    "\n",
+    "network_remapped = {}\n",
+    "for pattern, correction in network.items():\n",
+    "    pattern_mapped = tuple(flip_index_map_inverse[q] for q in pattern)\n",
+    "    correction_mapped = [correction_index_map_inverse[q] for q in correction]\n",
+    "    network_remapped[pattern_mapped] = correction_mapped\n",
+    "\n",
+    "network_remapped"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "python11",
+   "language": "python",
+   "name": "python11"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}

pyproject.toml

@@ -64,6 +64,8 @@ dependencies = [
   "qecsim>=1.0b9",
   "sinter>=1.14.0",
   "tqdm>=4.66.2",
+  "ortools>=9.14.2151",
+  "protobuf==6.31.1",
   "networkx>=3.4.2",
 ]
 dynamic = ["version"]
@@ -151,7 +153,7 @@ exclude = [
 
 [[tool.mypy.overrides]]
 module = ["qiskit.*", "qecsim.*", "qiskit_aer.*", "matplotlib.*", "scipy.*", "ldpc.*", "pytest_console_scripts.*",
-    "z3.*", "bposd.*", "numba.*", "pymatching.*", "stim.*", "multiprocess.*", "sinter.*", "qsample.*", "pandas.*", "tqdm.*", "networkx.*"]
+    "z3.*", "bposd.*", "numba.*", "pymatching.*", "stim.*", "multiprocess.*", "sinter.*", "qsample.*", "pandas.*", "tqdm.*", "ortools.*", "ortools.sat.python","networkx.*"]
 ignore_missing_imports = true
 
 
@@ -249,6 +251,7 @@ aer = "aer"
 anc = "anc"
 arange = "arange"
 nd = "nd"
+strat = "strat"
 iz = "iz"