Oxidize synth_permutation_acg (#12543)

* Move utility functions _inverse_pattern and _get_ordered_swap to Rust

* fix formatting and pylint issues

* Changed input type to `PyArrayLike1<i64, AllowTypeChange>`

* Refactor `permutation.rs`, clean up imports, fix coverage error

* fix docstring for `_inverse_pattern`

Co-authored-by: Raynel Sanchez <[email protected]>

* fix docstring for `_get_ordered_swap`

Co-authored-by: Raynel Sanchez <[email protected]>

* remove pymodule nesting

* remove explicit `AllowTypeChange`

* Move input validation out of `_inverse_pattern` and `_get_ordered_swap`

* oxidization attempt no. 1

* working version!

maybe faster possible...

* move more into rust & fix clones

* layouting & comments

* dangling comment

* circuit construction in rust

* remove dangling code

* more lint

* add reno

* drop redundant Ok(expect())

* Implements Shelly's suggestions

* simplify code a little

---------

Co-authored-by: jpacold <[email protected]>
Co-authored-by: Raynel Sanchez <[email protected]>

crates/accelerate/src/synthesis/permutation/mod.rs

@@ -59,10 +59,34 @@ pub fn _synth_permutation_basic(py: Python, pattern: PyArrayLike1<i64>) -> PyRes
     )
 }
 
+#[pyfunction]
+#[pyo3(signature = (pattern))]
+fn _synth_permutation_acg(py: Python, pattern: PyArrayLike1<i64>) -> PyResult<CircuitData> {
+    let inverted = utils::invert(&pattern.as_array());
+    let view = inverted.view();
+    let num_qubits = view.len();
+    let cycles = utils::pattern_to_cycles(&view);
+    let swaps = utils::decompose_cycles(&cycles);
+
+    CircuitData::from_standard_gates(
+        py,
+        num_qubits as u32,
+        swaps.iter().map(|(i, j)| {
+            (
+                StandardGate::SwapGate,
+                smallvec![],
+                smallvec![Qubit(*i as u32), Qubit(*j as u32)],
+            )
+        }),
+        Param::Float(0.0),
+    )
+}
+
 #[pymodule]
 pub fn permutation(m: &Bound<PyModule>) -> PyResult<()> {
     m.add_function(wrap_pyfunction!(_validate_permutation, m)?)?;
     m.add_function(wrap_pyfunction!(_inverse_pattern, m)?)?;
     m.add_function(wrap_pyfunction!(_synth_permutation_basic, m)?)?;
+    m.add_function(wrap_pyfunction!(_synth_permutation_acg, m)?)?;
     Ok(())
 }

crates/accelerate/src/synthesis/permutation/utils.rs

@@ -15,6 +15,8 @@ use pyo3::exceptions::PyValueError;
 use pyo3::prelude::*;
 use std::vec::Vec;
 
+use qiskit_circuit::slice::{PySequenceIndex, PySequenceIndexError, SequenceIndex};
+
 pub fn validate_permutation(pattern: &ArrayView1<i64>) -> PyResult<()> {
     let n = pattern.len();
     let mut seen: Vec<bool> = vec![false; n];
@@ -63,19 +65,19 @@ pub fn invert(pattern: &ArrayView1<i64>) -> Array1<usize> {
 /// then this creates a quantum circuit with ``m-1`` SWAPs (and of depth ``m-1``);
 /// if the input  permutation consists of several disjoint cycles, then each cycle
 /// is essentially treated independently.
-pub fn get_ordered_swap(pattern: &ArrayView1<i64>) -> Vec<(i64, i64)> {
+pub fn get_ordered_swap(pattern: &ArrayView1<i64>) -> Vec<(usize, usize)> {
     let mut permutation: Vec<usize> = pattern.iter().map(|&x| x as usize).collect();
     let mut index_map = invert(pattern);
 
     let n = permutation.len();
-    let mut swaps: Vec<(i64, i64)> = Vec::with_capacity(n);
+    let mut swaps: Vec<(usize, usize)> = Vec::with_capacity(n);
     for ii in 0..n {
         let val = permutation[ii];
         if val == ii {
             continue;
         }
         let jj = index_map[ii];
-        swaps.push((ii as i64, jj as i64));
+        swaps.push((ii, jj));
         (permutation[ii], permutation[jj]) = (permutation[jj], permutation[ii]);
         index_map[val] = jj;
         index_map[ii] = ii;
@@ -84,3 +86,68 @@ pub fn get_ordered_swap(pattern: &ArrayView1<i64>) -> Vec<(i64, i64)> {
     swaps[..].reverse();
     swaps
 }
+
+/// Explore cycles in a permutation pattern. This is probably best explained in an
+/// example: let a pattern be [1, 2, 3, 0, 4, 6, 5], then it contains the two
+/// cycles [1, 2, 3, 0] and [6, 5]. The index [4] does not perform a permutation and does
+/// therefore not create a cycle.
+pub fn pattern_to_cycles(pattern: &ArrayView1<usize>) -> Vec<Vec<usize>> {
+    // vector keeping track of which elements in the permutation pattern have been visited
+    let mut explored: Vec<bool> = vec![false; pattern.len()];
+
+    // vector to store the cycles
+    let mut cycles: Vec<Vec<usize>> = Vec::new();
+
+    for pos in pattern {
+        let mut cycle: Vec<usize> = Vec::new();
+
+        // follow the cycle until we reached an entry we saw before
+        let mut i = *pos;
+        while !explored[i] {
+            cycle.push(i);
+            explored[i] = true;
+            i = pattern[i];
+        }
+        // cycles must have more than 1 element
+        if cycle.len() > 1 {
+            cycles.push(cycle);
+        }
+    }
+
+    cycles
+}
+
+/// Periodic (or Python-like) access to a vector.
+/// Util used below in ``decompose_cycles``.
+#[inline]
+fn pget(vec: &Vec<usize>, index: isize) -> Result<usize, PySequenceIndexError> {
+    let SequenceIndex::Int(wrapped) = PySequenceIndex::Int(index).with_len(vec.len())? else {unreachable!()};
+    Ok(vec[wrapped])
+}
+
+/// Given a disjoint cycle decomposition of a permutation pattern (see the function
+/// ``pattern_to_cycles``), decomposes every cycle into a series of SWAPs to implement it.
+/// In combination with ``pattern_to_cycle``, this function allows to implement a
+/// full permutation pattern by applying SWAP gates on the returned index-pairs.
+pub fn decompose_cycles(cycles: &Vec<Vec<usize>>) -> Vec<(usize, usize)> {
+    let mut swaps: Vec<(usize, usize)> = Vec::new();
+
+    for cycle in cycles {
+        let length = cycle.len() as isize;
+
+        for idx in 0..(length - 1) / 2 {
+            swaps.push((
+                pget(cycle, idx - 1).unwrap(),
+                pget(cycle, length - 3 - idx).unwrap(),
+            ));
+        }
+        for idx in 0..length / 2 {
+            swaps.push((
+                pget(cycle, idx - 1).unwrap(),
+                pget(cycle, length - 2 - idx).unwrap(),
+            ));
+        }
+    }
+
+    swaps
+}

qiskit/synthesis/permutation/permutation_full.py

@@ -16,11 +16,9 @@
 
 import numpy as np
 from qiskit.circuit.quantumcircuit import QuantumCircuit
-from qiskit._accelerate.synthesis.permutation import _synth_permutation_basic
-from .permutation_utils import (
-    _inverse_pattern,
-    _pattern_to_cycles,
-    _decompose_cycles,
+from qiskit._accelerate.synthesis.permutation import (
+    _synth_permutation_basic,
+    _synth_permutation_acg,
 )
 
 
@@ -77,16 +75,4 @@ def synth_permutation_acg(pattern: list[int] | np.ndarray[int]) -> QuantumCircui
            *Routing Permutations on Graphs Via Matchings.*,
            `(Full paper) <https://www.cs.tau.ac.il/~nogaa/PDFS/r.pdf>`_
     """
-
-    num_qubits = len(pattern)
-    qc = QuantumCircuit(num_qubits)
-
-    # invert pattern (Qiskit notation is opposite)
-    cur_pattern = _inverse_pattern(pattern)
-    cycles = _pattern_to_cycles(cur_pattern)
-    swaps = _decompose_cycles(cycles)
-
-    for swap in swaps:
-        qc.swap(swap[0], swap[1])
-
-    return qc
+    return QuantumCircuit._from_circuit_data(_synth_permutation_acg(pattern))

qiskit/synthesis/permutation/permutation_utils.py

@@ -13,36 +13,4 @@
 """Utility functions for handling permutations."""
 
 # pylint: disable=unused-import
-from qiskit._accelerate.synthesis.permutation import (
-    _inverse_pattern,
-    _validate_permutation,
-)
-
-
-def _pattern_to_cycles(pattern):
-    """Given a permutation pattern, creates its disjoint cycle decomposition."""
-    nq = len(pattern)
-    explored = [False] * nq
-    cycles = []
-    for i in pattern:
-        cycle = []
-        while not explored[i]:
-            cycle.append(i)
-            explored[i] = True
-            i = pattern[i]
-        if len(cycle) >= 2:
-            cycles.append(cycle)
-    return cycles
-
-
-def _decompose_cycles(cycles):
-    """Given a disjoint cycle decomposition, decomposes every cycle into a SWAP
-    circuit of depth 2."""
-    swap_list = []
-    for cycle in cycles:
-        m = len(cycle)
-        for i in range((m - 1) // 2):
-            swap_list.append((cycle[i - 1], cycle[m - 3 - i]))
-        for i in range(m // 2):
-            swap_list.append((cycle[i - 1], cycle[m - 2 - i]))
-    return swap_list
+from qiskit._accelerate.synthesis.permutation import _inverse_pattern, _validate_permutation

releasenotes/notes/oxidize-acg-0294a87c0d5974fa.yaml

@@ -0,0 +1,5 @@
+---
+upgrade_synthesis:
+  - |
+    Port :func:`.synth_permutation_acg`, used to synthesize qubit permutations, to Rust.
+    This produces an approximate 3x performance improvement on 1000 qubit circuits.