Fully port Optimize1qGatesDecomposition to Rust (#12959)

* Fully port Optimize1qGatesDecomposition to Rust

This commit builds off of #12550 and the other data model in Rust
infrastructure and migrates the Optimize1qGatesDecomposition pass to
operate fully in Rust. The full path of the transpiler pass now never
leaves rust until it has finished modifying the DAGCircuit. There is
still some python interaction necessary to handle parts of the data
model that are still in Python, mainly calibrations and parameter
expressions (for global phase). But otherwise the entirety of the pass
operates in rust now.

This is just a first pass at the migration here, it moves the pass to be
a single for loop in rust. The next steps here are to look at operating
the pass in parallel. There is no data dependency between the
optimizations being done by the pass so we should be able to the
throughput of the pass by leveraging multithreading to handle each run
in parallel. This commit does not attempt this though, because of the
Python dependency and also the data structures around gates and the
dag aren't really setup for multithreading yet and there likely will
need to be some work to support that (this pass is a good candidate to
work through the bugs on that).

Part of #12208

* Tweak control_flow_op_nodes() method to avoid dag traversal when not necessary

* Store target basis set without heap allocation

Since we only are storing 12 enum fields (which are a single byte) using
any heap allocated collection is completely overkill and will have more
overhead that storing a statically sized array for all 12 variants. This
commit adds a new struct that wraps a `[bool; 12]` to track which
basis are supported and an API for tracking this. This simplifies the
tracking of which qubit supports which EulerBasis, it also means other
internal users of the 1q decomposition have a simplified API for working
with the euler basis.

* Remove From trait for Qubit->PhysicalQubit conversion

* Fix merge conflict

* Use new DAGCircuit::has_control_flow() for control_flow_op_nodes() pymethod

* Move _basis_gates set creation to __init__

* Update releasenotes/notes/optimize-1q-gates-decomposition-ce111961b6782ee0.yaml

Co-authored-by: Elena Peña Tapia <[email protected]>

---------

Co-authored-by: Elena Peña Tapia <[email protected]>

crates/accelerate/src/nlayout.rs

@@ -24,7 +24,7 @@ use hashbrown::HashMap;
 macro_rules! qubit_newtype {
     ($id: ident) => {
         #[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Hash)]
-        pub struct $id(u32);
+        pub struct $id(pub u32);
 
         impl $id {
             #[inline]
@@ -72,6 +72,7 @@ impl PhysicalQubit {
         layout.phys_to_virt[self.index()]
     }
 }
+
 qubit_newtype!(VirtualQubit);
 impl VirtualQubit {
     /// Get the physical qubit that currently corresponds to this index of virtual qubit in the

crates/accelerate/src/target_transpiler/mod.rs

@@ -940,6 +940,17 @@ impl Target {
         });
     }
 
+    /// Get the error rate of a given instruction in the target
+    pub fn get_error(&self, name: &str, qargs: &[PhysicalQubit]) -> Option<f64> {
+        self.gate_map.get(name).and_then(|gate_props| {
+            let qargs_key: Qargs = qargs.iter().cloned().collect();
+            match gate_props.get(Some(&qargs_key)) {
+                Some(props) => props.as_ref().and_then(|inst_props| inst_props.error),
+                None => None,
+            }
+        })
+    }
+
     /// Get an iterator over the indices of all physical qubits of the target
     pub fn physical_qubits(&self) -> impl ExactSizeIterator<Item = usize> {
         0..self.num_qubits.unwrap_or_default()

crates/accelerate/src/two_qubit_decompose.rs

@@ -41,7 +41,7 @@ use pyo3::types::PyList;
 
 use crate::convert_2q_block_matrix::change_basis;
 use crate::euler_one_qubit_decomposer::{
-    angles_from_unitary, det_one_qubit, unitary_to_gate_sequence_inner, EulerBasis,
+    angles_from_unitary, det_one_qubit, unitary_to_gate_sequence_inner, EulerBasis, EulerBasisSet,
     OneQubitGateSequence, ANGLE_ZERO_EPSILON,
 };
 use crate::utils;
@@ -1060,7 +1060,8 @@ impl TwoQubitWeylDecomposition {
             Some(basis) => EulerBasis::__new__(basis.deref())?,
             None => self.default_euler_basis,
         };
-        let target_1q_basis_list: Vec<EulerBasis> = vec![euler_basis];
+        let mut target_1q_basis_list = EulerBasisSet::new();
+        target_1q_basis_list.add_basis(euler_basis);
 
         let mut gate_sequence = CircuitData::with_capacity(py, 2, 0, 21, Param::Float(0.))?;
         let mut global_phase: f64 = self.global_phase;
@@ -1507,7 +1508,8 @@ impl TwoQubitBasisDecomposer {
         unitary: ArrayView2<Complex64>,
         qubit: u8,
     ) {
-        let target_1q_basis_list = vec![self.euler_basis];
+        let mut target_1q_basis_list = EulerBasisSet::new();
+        target_1q_basis_list.add_basis(self.euler_basis);
         let sequence = unitary_to_gate_sequence_inner(
             unitary,
             &target_1q_basis_list,
@@ -1751,7 +1753,8 @@ impl TwoQubitBasisDecomposer {
         if let Some(seq) = sequence {
             return Ok(seq);
         }
-        let target_1q_basis_list = vec![self.euler_basis];
+        let mut target_1q_basis_list = EulerBasisSet::new();
+        target_1q_basis_list.add_basis(self.euler_basis);
         let euler_decompositions: SmallVec<[Option<OneQubitGateSequence>; 8]> = decomposition
             .iter()
             .map(|decomp| {
@@ -1993,7 +1996,8 @@ impl TwoQubitBasisDecomposer {
         if let Some(seq) = sequence {
             return Ok(seq);
         }
-        let target_1q_basis_list = vec![self.euler_basis];
+        let mut target_1q_basis_list = EulerBasisSet::new();
+        target_1q_basis_list.add_basis(self.euler_basis);
         let euler_decompositions: SmallVec<[Option<OneQubitGateSequence>; 8]> = decomposition
             .iter()
             .map(|decomp| {

crates/circuit/src/dag_circuit.rs

@@ -22,8 +22,8 @@ use crate::dag_node::{DAGInNode, DAGNode, DAGOpNode, DAGOutNode};
 use crate::dot_utils::build_dot;
 use crate::error::DAGCircuitError;
 use crate::imports;
-use crate::interner::Interner;
-use crate::operations::{Operation, OperationRef, Param, PyInstruction};
+use crate::interner::{Interned, Interner};
+use crate::operations::{Operation, OperationRef, Param, PyInstruction, StandardGate};
 use crate::packed_instruction::PackedInstruction;
 use crate::rustworkx_core_vnext::isomorphism;
 use crate::{BitType, Clbit, Qubit, TupleLikeArg};
@@ -3867,7 +3867,7 @@ def _format(operand):
     ///     include_directives (bool): include `barrier`, `snapshot` etc.
     ///
     /// Returns:
-    ///     list[DAGOpNode]: the list of node ids containing the given op.
+    ///     list[DAGOpNode]: the list of dag nodes containing the given op.
     #[pyo3(name= "op_nodes", signature=(op=None, include_directives=true))]
     fn py_op_nodes(
         &self,
@@ -3903,6 +3903,33 @@ def _format(operand):
         Ok(nodes)
     }
 
+    /// Get a list of "op" nodes in the dag that contain control flow instructions.
+    ///
+    /// Returns:
+    ///     list[DAGOpNode] | None: The list of dag nodes containing control flow ops. If there
+    ///         are no control flow nodes None is returned
+    fn control_flow_op_nodes(&self, py: Python) -> PyResult<Option<Vec<Py<PyAny>>>> {
+        if self.has_control_flow() {
+            let result: PyResult<Vec<Py<PyAny>>> = self
+                .dag
+                .node_references()
+                .filter_map(|(node_index, node_type)| match node_type {
+                    NodeType::Operation(ref node) => {
+                        if node.op.control_flow() {
+                            Some(self.unpack_into(py, node_index, node_type))
+                        } else {
+                            None
+                        }
+                    }
+                    _ => None,
+                })
+                .collect();
+            Ok(Some(result?))
+        } else {
+            Ok(None)
+        }
+    }
+
     /// Get the list of gate nodes in the dag.
     ///
     /// Returns:
@@ -4978,31 +5005,6 @@ def _format(operand):
         Ok(result)
     }
 
-    fn _insert_1q_on_incoming_qubit(
-        &mut self,
-        py: Python,
-        node: &Bound<PyAny>,
-        old_index: usize,
-    ) -> PyResult<()> {
-        if let NodeType::Operation(inst) = self.pack_into(py, node)? {
-            self.increment_op(inst.op.name());
-            let new_index = self.dag.add_node(NodeType::Operation(inst));
-            let old_index: NodeIndex = NodeIndex::new(old_index);
-            let (parent_index, edge_index, weight) = self
-                .dag
-                .edges_directed(old_index, Incoming)
-                .map(|edge| (edge.source(), edge.id(), edge.weight().clone()))
-                .next()
-                .unwrap();
-            self.dag.add_edge(parent_index, new_index, weight.clone());
-            self.dag.add_edge(new_index, old_index, weight);
-            self.dag.remove_edge(edge_index);
-            Ok(())
-        } else {
-            Err(PyTypeError::new_err("Invalid node type input"))
-        }
-    }
-
     fn _edges(&self, py: Python) -> Vec<PyObject> {
         self.dag
             .edge_indices()
@@ -5604,7 +5606,7 @@ impl DAGCircuit {
     /// Remove an operation node n.
     ///
     /// Add edges from predecessors to successors.
-    fn remove_op_node(&mut self, index: NodeIndex) {
+    pub fn remove_op_node(&mut self, index: NodeIndex) {
         let mut edge_list: Vec<(NodeIndex, NodeIndex, Wire)> = Vec::new();
         for (source, in_weight) in self
             .dag
@@ -6154,6 +6156,127 @@ impl DAGCircuit {
         }
         Ok(())
     }
+    /// Get qargs from an intern index
+    pub fn get_qargs(&self, index: Interned<[Qubit]>) -> &[Qubit] {
+        self.qargs_interner.get(index)
+    }
+
+    /// Get cargs from an intern index
+    pub fn get_cargs(&self, index: Interned<[Clbit]>) -> &[Clbit] {
+        self.cargs_interner.get(index)
+    }
+
+    /// Insert a new 1q standard gate on incoming qubit
+    pub fn insert_1q_on_incoming_qubit(
+        &mut self,
+        new_gate: (StandardGate, &[f64]),
+        old_index: NodeIndex,
+    ) {
+        self.increment_op(new_gate.0.name());
+        let old_node = &self.dag[old_index];
+        let inst = if let NodeType::Operation(old_node) = old_node {
+            PackedInstruction {
+                op: new_gate.0.into(),
+                qubits: old_node.qubits,
+                clbits: old_node.clbits,
+                params: (!new_gate.1.is_empty())
+                    .then(|| Box::new(new_gate.1.iter().map(|x| Param::Float(*x)).collect())),
+                extra_attrs: None,
+                #[cfg(feature = "cache_pygates")]
+                py_op: OnceCell::new(),
+            }
+        } else {
+            panic!("This method only works if provided index is an op node");
+        };
+        let new_index = self.dag.add_node(NodeType::Operation(inst));
+        let (parent_index, edge_index, weight) = self
+            .dag
+            .edges_directed(old_index, Incoming)
+            .map(|edge| (edge.source(), edge.id(), edge.weight().clone()))
+            .next()
+            .unwrap();
+        self.dag.add_edge(parent_index, new_index, weight.clone());
+        self.dag.add_edge(new_index, old_index, weight);
+        self.dag.remove_edge(edge_index);
+    }
+
+    /// Remove a sequence of 1 qubit nodes from the dag
+    /// This must only be called if all the nodes operate
+    /// on a single qubit with no other wires in or out of any nodes
+    pub fn remove_1q_sequence(&mut self, sequence: &[NodeIndex]) {
+        let (parent_index, weight) = self
+            .dag
+            .edges_directed(*sequence.first().unwrap(), Incoming)
+            .map(|edge| (edge.source(), edge.weight().clone()))
+            .next()
+            .unwrap();
+        let child_index = self
+            .dag
+            .edges_directed(*sequence.last().unwrap(), Outgoing)
+            .map(|edge| edge.target())
+            .next()
+            .unwrap();
+        self.dag.add_edge(parent_index, child_index, weight);
+        for node in sequence {
+            match self.dag.remove_node(*node) {
+                Some(NodeType::Operation(packed)) => {
+                    let op_name = packed.op.name();
+                    self.decrement_op(op_name);
+                }
+                _ => panic!("Must be called with valid operation node!"),
+            }
+        }
+    }
+
+    pub fn add_global_phase(&mut self, py: Python, value: &Param) -> PyResult<()> {
+        match value {
+            Param::Obj(_) => {
+                return Err(PyTypeError::new_err(
+                    "Invalid parameter type, only float and parameter expression are supported",
+                ))
+            }
+            _ => self.set_global_phase(add_global_phase(py, &self.global_phase, value)?)?,
+        }
+        Ok(())
+    }
+
+    pub fn calibrations_empty(&self) -> bool {
+        self.calibrations.is_empty()
+    }
+
+    pub fn has_calibration_for_index(&self, py: Python, node_index: NodeIndex) -> PyResult<bool> {
+        let node = &self.dag[node_index];
+        if let NodeType::Operation(instruction) = node {
+            if !self.calibrations.contains_key(instruction.op.name()) {
+                return Ok(false);
+            }
+            let params = match &instruction.params {
+                Some(params) => {
+                    let mut out_params = Vec::new();
+                    for p in params.iter() {
+                        if let Param::ParameterExpression(exp) = p {
+                            let exp = exp.bind(py);
+                            if !exp.getattr(intern!(py, "parameters"))?.is_truthy()? {
+                                let as_py_float = exp.call_method0(intern!(py, "__float__"))?;
+                                out_params.push(as_py_float.unbind());
+                                continue;
+                            }
+                        }
+                        out_params.push(p.to_object(py));
+                    }
+                    PyTuple::new_bound(py, out_params)
+                }
+                None => PyTuple::empty_bound(py),
+            };
+            let qargs = self.qargs_interner.get(instruction.qubits);
+            let qubits = PyTuple::new_bound(py, qargs.iter().map(|x| x.0));
+            self.calibrations[instruction.op.name()]
+                .bind(py)
+                .contains((qubits, params).to_object(py))
+        } else {
+            Err(DAGCircuitError::new_err("Specified node is not an op node"))
+        }
+    }
 }
 
 /// Add to global phase. Global phase can only be Float or ParameterExpression so this

qiskit/transpiler/passes/optimization/optimize_1q_decomposition.py

@@ -35,7 +35,6 @@
 from qiskit.circuit import Qubit
 from qiskit.circuit.quantumcircuitdata import CircuitInstruction
 from qiskit.dagcircuit.dagcircuit import DAGCircuit
-from qiskit.dagcircuit.dagnode import DAGOpNode
 
 
 logger = logging.getLogger(__name__)
@@ -82,9 +81,12 @@ def __init__(self, basis=None, target=None):
         """
         super().__init__()
 
-        self._basis_gates = basis
+        if basis:
+            self._basis_gates = set(basis)
+        else:
+            self._basis_gates = None
         self._target = target
-        self._global_decomposers = []
+        self._global_decomposers = None
         self._local_decomposers_cache = {}
 
         if basis:
@@ -209,46 +211,12 @@ def run(self, dag):
         Returns:
             DAGCircuit: the optimized DAG.
         """
-        runs = []
-        qubits = []
-        bases = []
-        for run in dag.collect_1q_runs():
-            qubit = dag.find_bit(run[0].qargs[0]).index
-            runs.append(run)
-            qubits.append(qubit)
-            bases.append(self._get_decomposer(qubit))
-        best_sequences = euler_one_qubit_decomposer.optimize_1q_gates_decomposition(
-            runs, qubits, bases, simplify=True, error_map=self.error_map
+        euler_one_qubit_decomposer.optimize_1q_gates_decomposition(
+            dag,
+            target=self._target,
+            global_decomposers=self._global_decomposers,
+            basis_gates=self._basis_gates,
         )
-        for index, best_circuit_sequence in enumerate(best_sequences):
-            run = runs[index]
-            qubit = qubits[index]
-            if self._target is None:
-                basis = self._basis_gates
-            else:
-                basis = self._target.operation_names_for_qargs((qubit,))
-            if best_circuit_sequence is not None:
-                (old_error, new_error, best_circuit_sequence) = best_circuit_sequence
-                if self._substitution_checks(
-                    dag,
-                    run,
-                    best_circuit_sequence,
-                    basis,
-                    qubit,
-                    old_error=old_error,
-                    new_error=new_error,
-                ):
-                    first_node_id = run[0]._node_id
-                    qubit = run[0].qargs
-                    for gate, angles in best_circuit_sequence:
-                        op = CircuitInstruction.from_standard(gate, qubit, angles)
-                        node = DAGOpNode.from_instruction(op)
-                        dag._insert_1q_on_incoming_qubit(node, first_node_id)
-                    dag.global_phase += best_circuit_sequence.global_phase
-                    # Delete the other nodes in the run
-                    for current_node in run:
-                        dag.remove_op_node(current_node)
-
         return dag
 
     def _error(self, circuit, qubit):

qiskit/transpiler/passes/utils/control_flow.py

@@ -54,10 +54,12 @@ def out(self, dag):
         def bound_wrapped_method(dag):
             return out(self, dag)
 
-        for node in dag.op_nodes(ControlFlowOp):
-            dag.substitute_node(
-                node, map_blocks(bound_wrapped_method, node.op), propagate_condition=False
-            )
+        control_flow_nodes = dag.control_flow_op_nodes()
+        if control_flow_nodes is not None:
+            for node in control_flow_nodes:
+                dag.substitute_node(
+                    node, map_blocks(bound_wrapped_method, node.op), propagate_condition=False
+                )
         return method(self, dag)
 
     return out