From bfc69a499c035bc3b8af31b37ba0469d7e89e67a Mon Sep 17 00:00:00 2001 From: Matthew Treinish Date: Wed, 3 Jul 2024 00:18:40 -0400 Subject: [PATCH] Use rust gates for Optimize1QGatesDecomposition (#12650) * Use rust gates for Optimize1QGatesDecomposition This commit moves to using rust gates for the Optimize1QGatesDecomposition transpiler pass. It takes in a sequence of runs (which are a list of DAGOpNodes) from the python side of the transpiler pass which are generated from DAGCircuit.collect_1q_runs() (which in the future should be moved to rust after #12550 merges). The rust portion of the pass now iterates over each run, performs the matrix multiplication to compute the unitary of the run, then synthesizes that unitary, computes the estimated error of the circuit synthesis and returns a tuple of the circuit sequence in terms of rust StandardGate enums. The python portion of the code then takes those sequences and does inplace substitution of each run with the sequence returned from rust. Once #12550 merges we should be able to move the input collect_1q_runs() call and perform the output node substitions in rust making the full pass execute in the rust domain without any python interaction. Additionally, the OneQubitEulerDecomposer class is updated to use rust for circuit generation instead of doing this python side. The internal changes done to use rust gates in the transpiler pass meant we were half way to this already by emitting rust StandardGates instead of python gate objects. The dag handling is still done in Python however until #12550 merges. This also includes an implementation of the r gate, I temporarily added this to unblock this effort as it was the only gate missing needed to complete this. We can rebase this if a standalone implementation of the gate merges before this. * Cache target decompositions for each qubit Previously this PR was re-computing the target bases to synthesize with for each run found in the circuit. But in cases where there were multiple runs repeated on a qubit this was unecessary work. Prior to moving this code to rust there was already caching code to make this optimization, but the rust path short circuited around this. This commit fixes this so we're caching the target bases for each qubit and only computing it once. * Optimize rust implementation slightly * Avoid extra allocations by inlining matrix multiplication * Remove unnecessary comment * Remove stray code block * Add import path for rust gate * Use rust gate in circuit constructor * Remove duplicated op_name getter and just use existing name getter * Apply suggestions from code review Co-authored-by: John Lapeyre * Simplify construction of target_basis_vec * Fix rebase issue * Update crates/accelerate/src/euler_one_qubit_decomposer.rs Co-authored-by: John Lapeyre * Update crates/accelerate/src/euler_one_qubit_decomposer.rs --------- Co-authored-by: John Lapeyre --- .../src/euler_one_qubit_decomposer.rs | 301 +++++++++++++++--- crates/accelerate/src/two_qubit_decompose.rs | 19 +- crates/circuit/src/dag_node.rs | 50 ++- crates/circuit/src/operations.rs | 6 + qiskit/dagcircuit/dagcircuit.py | 57 ++-- .../one_qubit/one_qubit_decompose.py | 50 ++- .../optimization/optimize_1q_decomposition.py | 105 +++--- 7 files changed, 431 insertions(+), 157 deletions(-) diff --git a/crates/accelerate/src/euler_one_qubit_decomposer.rs b/crates/accelerate/src/euler_one_qubit_decomposer.rs index 01725269bb84..8c3a87ce51ec 100644 --- a/crates/accelerate/src/euler_one_qubit_decomposer.rs +++ b/crates/accelerate/src/euler_one_qubit_decomposer.rs @@ -18,7 +18,6 @@ use num_complex::{Complex64, ComplexFloat}; use smallvec::{smallvec, SmallVec}; use std::cmp::Ordering; use std::f64::consts::PI; -use std::ops::Deref; use std::str::FromStr; use pyo3::exceptions::PyValueError; @@ -31,8 +30,12 @@ use ndarray::prelude::*; use numpy::PyReadonlyArray2; use pyo3::pybacked::PyBackedStr; +use qiskit_circuit::circuit_data::CircuitData; +use qiskit_circuit::dag_node::DAGOpNode; +use qiskit_circuit::operations::{Operation, Param, StandardGate}; use qiskit_circuit::slice::{PySequenceIndex, SequenceIndex}; use qiskit_circuit::util::c64; +use qiskit_circuit::Qubit; pub const ANGLE_ZERO_EPSILON: f64 = 1e-12; @@ -68,12 +71,12 @@ impl OneQubitGateErrorMap { #[pyclass(sequence)] pub struct OneQubitGateSequence { - pub gates: Vec<(String, SmallVec<[f64; 3]>)>, + pub gates: Vec<(StandardGate, SmallVec<[f64; 3]>)>, #[pyo3(get)] pub global_phase: f64, } -type OneQubitGateSequenceState = (Vec<(String, SmallVec<[f64; 3]>)>, f64); +type OneQubitGateSequenceState = (Vec<(StandardGate, SmallVec<[f64; 3]>)>, f64); #[pymethods] impl OneQubitGateSequence { @@ -115,15 +118,15 @@ fn circuit_kak( phi: f64, lam: f64, phase: f64, - k_gate: &str, - a_gate: &str, + k_gate: StandardGate, + a_gate: StandardGate, simplify: bool, atol: Option, ) -> OneQubitGateSequence { let mut lam = lam; let mut theta = theta; let mut phi = phi; - let mut circuit: Vec<(String, SmallVec<[f64; 3]>)> = Vec::with_capacity(3); + let mut circuit: Vec<(StandardGate, SmallVec<[f64; 3]>)> = Vec::with_capacity(3); let mut atol = match atol { Some(atol) => atol, None => ANGLE_ZERO_EPSILON, @@ -139,7 +142,7 @@ fn circuit_kak( // slippage coming from _mod_2pi injecting multiples of 2pi. lam = mod_2pi(lam, atol); if lam.abs() > atol { - circuit.push((String::from(k_gate), smallvec![lam])); + circuit.push((k_gate, smallvec![lam])); global_phase += lam / 2.; } return OneQubitGateSequence { @@ -160,13 +163,13 @@ fn circuit_kak( lam = mod_2pi(lam, atol); if lam.abs() > atol { global_phase += lam / 2.; - circuit.push((String::from(k_gate), smallvec![lam])); + circuit.push((k_gate, smallvec![lam])); } - circuit.push((String::from(a_gate), smallvec![theta])); + circuit.push((a_gate, smallvec![theta])); phi = mod_2pi(phi, atol); if phi.abs() > atol { global_phase += phi / 2.; - circuit.push((String::from(k_gate), smallvec![phi])); + circuit.push((k_gate, smallvec![phi])); } OneQubitGateSequence { gates: circuit, @@ -190,7 +193,7 @@ fn circuit_u3( let phi = mod_2pi(phi, atol); let lam = mod_2pi(lam, atol); if !simplify || theta.abs() > atol || phi.abs() > atol || lam.abs() > atol { - circuit.push((String::from("u3"), smallvec![theta, phi, lam])); + circuit.push((StandardGate::U3Gate, smallvec![theta, phi, lam])); } OneQubitGateSequence { gates: circuit, @@ -217,16 +220,16 @@ fn circuit_u321( if theta.abs() < atol { let tot = mod_2pi(phi + lam, atol); if tot.abs() > atol { - circuit.push((String::from("u1"), smallvec![tot])); + circuit.push((StandardGate::U1Gate, smallvec![tot])); } } else if (theta - PI / 2.).abs() < atol { circuit.push(( - String::from("u2"), + StandardGate::U2Gate, smallvec![mod_2pi(phi, atol), mod_2pi(lam, atol)], )); } else { circuit.push(( - String::from("u3"), + StandardGate::U3Gate, smallvec![theta, mod_2pi(phi, atol), mod_2pi(lam, atol)], )); } @@ -255,7 +258,7 @@ fn circuit_u( let phi = mod_2pi(phi, atol); let lam = mod_2pi(lam, atol); if theta.abs() > atol || phi.abs() > atol || lam.abs() > atol { - circuit.push((String::from("u"), smallvec![theta, phi, lam])); + circuit.push((StandardGate::UGate, smallvec![theta, phi, lam])); } OneQubitGateSequence { gates: circuit, @@ -358,7 +361,7 @@ fn circuit_rr( // This can be expressed as a single R gate if theta.abs() > atol { circuit.push(( - String::from("r"), + StandardGate::RGate, smallvec![theta, mod_2pi(PI / 2. + phi, atol)], )); } @@ -366,12 +369,12 @@ fn circuit_rr( // General case: use two R gates if (theta - PI).abs() > atol { circuit.push(( - String::from("r"), + StandardGate::RGate, smallvec![theta - PI, mod_2pi(PI / 2. - lam, atol)], )); } circuit.push(( - String::from("r"), + StandardGate::RGate, smallvec![PI, mod_2pi(0.5 * (phi - lam + PI), atol)], )); } @@ -393,10 +396,46 @@ pub fn generate_circuit( atol: Option, ) -> PyResult { let res = match target_basis { - EulerBasis::ZYZ => circuit_kak(theta, phi, lam, phase, "rz", "ry", simplify, atol), - EulerBasis::ZXZ => circuit_kak(theta, phi, lam, phase, "rz", "rx", simplify, atol), - EulerBasis::XZX => circuit_kak(theta, phi, lam, phase, "rx", "rz", simplify, atol), - EulerBasis::XYX => circuit_kak(theta, phi, lam, phase, "rx", "ry", simplify, atol), + EulerBasis::ZYZ => circuit_kak( + theta, + phi, + lam, + phase, + StandardGate::RZGate, + StandardGate::RYGate, + simplify, + atol, + ), + EulerBasis::ZXZ => circuit_kak( + theta, + phi, + lam, + phase, + StandardGate::RZGate, + StandardGate::RXGate, + simplify, + atol, + ), + EulerBasis::XZX => circuit_kak( + theta, + phi, + lam, + phase, + StandardGate::RXGate, + StandardGate::RZGate, + simplify, + atol, + ), + EulerBasis::XYX => circuit_kak( + theta, + phi, + lam, + phase, + StandardGate::RXGate, + StandardGate::RYGate, + simplify, + atol, + ), EulerBasis::U3 => circuit_u3(theta, phi, lam, phase, simplify, atol), EulerBasis::U321 => circuit_u321(theta, phi, lam, phase, simplify, atol), EulerBasis::U => circuit_u(theta, phi, lam, phase, simplify, atol), @@ -411,11 +450,13 @@ pub fn generate_circuit( let fnz = |circuit: &mut OneQubitGateSequence, phi: f64| { let phi = mod_2pi(phi, inner_atol); if phi.abs() > inner_atol { - circuit.gates.push((String::from("p"), smallvec![phi])); + circuit + .gates + .push((StandardGate::PhaseGate, smallvec![phi])); } }; let fnx = |circuit: &mut OneQubitGateSequence| { - circuit.gates.push((String::from("sx"), SmallVec::new())); + circuit.gates.push((StandardGate::SXGate, SmallVec::new())); }; circuit_psx_gen( @@ -441,12 +482,12 @@ pub fn generate_circuit( let fnz = |circuit: &mut OneQubitGateSequence, phi: f64| { let phi = mod_2pi(phi, inner_atol); if phi.abs() > inner_atol { - circuit.gates.push((String::from("rz"), smallvec![phi])); + circuit.gates.push((StandardGate::RZGate, smallvec![phi])); circuit.global_phase += phi / 2.; } }; let fnx = |circuit: &mut OneQubitGateSequence| { - circuit.gates.push((String::from("sx"), SmallVec::new())); + circuit.gates.push((StandardGate::SXGate, SmallVec::new())); }; circuit_psx_gen( theta, @@ -471,12 +512,14 @@ pub fn generate_circuit( let fnz = |circuit: &mut OneQubitGateSequence, phi: f64| { let phi = mod_2pi(phi, inner_atol); if phi.abs() > inner_atol { - circuit.gates.push((String::from("u1"), smallvec![phi])); + circuit.gates.push((StandardGate::U1Gate, smallvec![phi])); } }; let fnx = |circuit: &mut OneQubitGateSequence| { circuit.global_phase += PI / 4.; - circuit.gates.push((String::from("rx"), smallvec![PI / 2.])); + circuit + .gates + .push((StandardGate::RXGate, smallvec![PI / 2.])); }; circuit_psx_gen( theta, @@ -501,15 +544,15 @@ pub fn generate_circuit( let fnz = |circuit: &mut OneQubitGateSequence, phi: f64| { let phi = mod_2pi(phi, inner_atol); if phi.abs() > inner_atol { - circuit.gates.push((String::from("rz"), smallvec![phi])); + circuit.gates.push((StandardGate::RZGate, smallvec![phi])); circuit.global_phase += phi / 2.; } }; let fnx = |circuit: &mut OneQubitGateSequence| { - circuit.gates.push((String::from("sx"), SmallVec::new())); + circuit.gates.push((StandardGate::SXGate, SmallVec::new())); }; let fnxpi = |circuit: &mut OneQubitGateSequence| { - circuit.gates.push((String::from("x"), SmallVec::new())); + circuit.gates.push((StandardGate::XGate, SmallVec::new())); }; circuit_psx_gen( theta, @@ -633,7 +676,7 @@ fn compare_error_fn( let fidelity_product: f64 = circuit .gates .iter() - .map(|x| 1. - err_map.get(&x.0).unwrap_or(&0.)) + .map(|gate| 1. - err_map.get(gate.0.name()).unwrap_or(&0.)) .product(); (1. - fidelity_product, circuit.gates.len()) } @@ -642,6 +685,28 @@ fn compare_error_fn( } fn compute_error( + gates: &[(StandardGate, SmallVec<[f64; 3]>)], + error_map: Option<&OneQubitGateErrorMap>, + qubit: usize, +) -> (f64, usize) { + match error_map { + Some(err_map) => { + let num_gates = gates.len(); + let gate_fidelities: f64 = gates + .iter() + .map(|gate| 1. - err_map.error_map[qubit].get(gate.0.name()).unwrap_or(&0.)) + .product(); + (1. - gate_fidelities, num_gates) + } + None => (gates.len() as f64, gates.len()), + } +} + +fn compute_error_term(gate: &str, error_map: &OneQubitGateErrorMap, qubit: usize) -> f64 { + 1. - error_map.error_map[qubit].get(gate).unwrap_or(&0.) +} + +fn compute_error_str( gates: &[(String, SmallVec<[f64; 3]>)], error_map: Option<&OneQubitGateErrorMap>, qubit: usize, @@ -651,7 +716,7 @@ fn compute_error( let num_gates = gates.len(); let gate_fidelities: f64 = gates .iter() - .map(|x| 1. - err_map.error_map[qubit].get(&x.0).unwrap_or(&0.)) + .map(|gate| compute_error_term(gate.0.as_str(), err_map, qubit)) .product(); (1. - gate_fidelities, num_gates) } @@ -670,11 +735,20 @@ pub fn compute_error_one_qubit_sequence( #[pyfunction] pub fn compute_error_list( - circuit: Vec<(String, SmallVec<[f64; 3]>)>, + circuit: Vec>, qubit: usize, error_map: Option<&OneQubitGateErrorMap>, ) -> (f64, usize) { - compute_error(&circuit, error_map, qubit) + let circuit_list: Vec<(String, SmallVec<[f64; 3]>)> = circuit + .iter() + .map(|node| { + ( + node.instruction.operation.name().to_string(), + smallvec![], // Params not needed in this path + ) + }) + .collect(); + compute_error_str(&circuit_list, error_map, qubit) } #[pyfunction] @@ -687,15 +761,13 @@ pub fn unitary_to_gate_sequence( simplify: bool, atol: Option, ) -> PyResult> { - let mut target_basis_vec: Vec = Vec::with_capacity(target_basis_list.len()); - for basis in target_basis_list { - let basis_enum = EulerBasis::__new__(basis.deref())?; - target_basis_vec.push(basis_enum) - } - let unitary_mat = unitary.as_array(); + let target_basis_vec: PyResult> = target_basis_list + .iter() + .map(|basis| EulerBasis::__new__(basis)) + .collect(); Ok(unitary_to_gate_sequence_inner( - unitary_mat, - &target_basis_vec, + unitary.as_array(), + &target_basis_vec?, qubit, error_map, simplify, @@ -725,6 +797,46 @@ pub fn unitary_to_gate_sequence_inner( }) } +#[pyfunction] +#[pyo3(signature = (unitary, target_basis_list, qubit, error_map=None, simplify=true, atol=None))] +pub fn unitary_to_circuit( + py: Python, + unitary: PyReadonlyArray2, + target_basis_list: Vec, + qubit: usize, + error_map: Option<&OneQubitGateErrorMap>, + simplify: bool, + atol: Option, +) -> PyResult> { + let target_basis_vec: PyResult> = target_basis_list + .iter() + .map(|basis| EulerBasis::__new__(basis)) + .collect(); + let circuit_sequence = unitary_to_gate_sequence_inner( + unitary.as_array(), + &target_basis_vec?, + qubit, + error_map, + simplify, + atol, + ); + Ok(circuit_sequence.map(|seq| { + CircuitData::from_standard_gates( + py, + 1, + seq.gates.into_iter().map(|(gate, params)| { + ( + gate, + params.into_iter().map(Param::Float).collect(), + smallvec![Qubit(0)], + ) + }), + Param::Float(seq.global_phase), + ) + .expect("Unexpected Qiskit python bug") + })) +} + #[inline] pub fn det_one_qubit(mat: ArrayView2) -> Complex64 { mat[[0, 0]] * mat[[1, 1]] - mat[[0, 1]] * mat[[1, 0]] @@ -853,6 +965,106 @@ pub fn params_zxz(unitary: PyReadonlyArray2) -> [f64; 4] { params_zxz_inner(mat) } +type OptimizeDecompositionReturn = Option<((f64, usize), (f64, usize), OneQubitGateSequence)>; + +#[pyfunction] +pub fn optimize_1q_gates_decomposition( + runs: Vec>>, + qubits: Vec, + bases: Vec>, + simplify: bool, + error_map: Option<&OneQubitGateErrorMap>, + atol: Option, +) -> Vec { + runs.iter() + .enumerate() + .map(|(index, raw_run)| -> OptimizeDecompositionReturn { + let mut error = match error_map { + Some(_) => 1., + None => raw_run.len() as f64, + }; + let qubit = qubits[index]; + let operator = &raw_run + .iter() + .map(|node| { + if let Some(err_map) = error_map { + error *= + compute_error_term(node.instruction.operation.name(), err_map, qubit) + } + node.instruction + .operation + .matrix(&node.instruction.params) + .expect("No matrix defined for operation") + }) + .fold( + [ + [Complex64::new(1., 0.), Complex64::new(0., 0.)], + [Complex64::new(0., 0.), Complex64::new(1., 0.)], + ], + |mut operator, node| { + matmul_1q(&mut operator, node); + operator + }, + ); + let old_error = if error_map.is_some() { + (1. - error, raw_run.len()) + } else { + (error, raw_run.len()) + }; + let target_basis_vec: Vec = bases[index] + .iter() + .map(|basis| EulerBasis::__new__(basis).unwrap()) + .collect(); + unitary_to_gate_sequence_inner( + aview2(operator), + &target_basis_vec, + qubit, + error_map, + simplify, + atol, + ) + .map(|out_seq| { + let new_error = compute_error_one_qubit_sequence(&out_seq, qubit, error_map); + (old_error, new_error, out_seq) + }) + }) + .collect() +} + +fn matmul_1q(operator: &mut [[Complex64; 2]; 2], other: Array2) { + *operator = [ + [ + other[[0, 0]] * operator[0][0] + other[[0, 1]] * operator[1][0], + other[[0, 0]] * operator[0][1] + other[[0, 1]] * operator[1][1], + ], + [ + other[[1, 0]] * operator[0][0] + other[[1, 1]] * operator[1][0], + other[[1, 0]] * operator[0][1] + other[[1, 1]] * operator[1][1], + ], + ]; +} + +#[pyfunction] +pub fn collect_1q_runs_filter(py: Python, node: PyObject) -> bool { + let op_node = node.extract::>(py); + match op_node { + Ok(node) => { + node.instruction.operation.num_qubits() == 1 + && node.instruction.operation.num_clbits() == 0 + && node + .instruction + .operation + .matrix(&node.instruction.params) + .is_some() + && match &node.instruction.extra_attrs { + None => true, + Some(attrs) => attrs.condition.is_none(), + } + } + Err(_) => false, + } +} + #[pymodule] pub fn euler_one_qubit_decomposer(m: &Bound) -> PyResult<()> { m.add_wrapped(wrap_pyfunction!(params_zyz))?; @@ -863,8 +1075,11 @@ pub fn euler_one_qubit_decomposer(m: &Bound) -> PyResult<()> { m.add_wrapped(wrap_pyfunction!(params_u1x))?; m.add_wrapped(wrap_pyfunction!(generate_circuit))?; m.add_wrapped(wrap_pyfunction!(unitary_to_gate_sequence))?; + m.add_wrapped(wrap_pyfunction!(unitary_to_circuit))?; m.add_wrapped(wrap_pyfunction!(compute_error_one_qubit_sequence))?; m.add_wrapped(wrap_pyfunction!(compute_error_list))?; + m.add_wrapped(wrap_pyfunction!(optimize_1q_gates_decomposition))?; + m.add_wrapped(wrap_pyfunction!(collect_1q_runs_filter))?; m.add_class::()?; m.add_class::()?; m.add_class::()?; diff --git a/crates/accelerate/src/two_qubit_decompose.rs b/crates/accelerate/src/two_qubit_decompose.rs index 37061d5159f4..568206925c2b 100644 --- a/crates/accelerate/src/two_qubit_decompose.rs +++ b/crates/accelerate/src/two_qubit_decompose.rs @@ -52,6 +52,7 @@ use rand_distr::StandardNormal; use rand_pcg::Pcg64Mcg; use qiskit_circuit::gate_matrix::{CX_GATE, H_GATE, ONE_QUBIT_IDENTITY, SX_GATE, X_GATE}; +use qiskit_circuit::operations::Operation; use qiskit_circuit::slice::{PySequenceIndex, SequenceIndex}; use qiskit_circuit::util::{c64, GateArray1Q, GateArray2Q, C_M_ONE, C_ONE, C_ZERO, IM, M_IM}; @@ -1045,7 +1046,7 @@ impl TwoQubitWeylDecomposition { ) .unwrap(); for gate in c2r.gates { - gate_sequence.push((gate.0, gate.1, smallvec![0])) + gate_sequence.push((gate.0.name().to_string(), gate.1, smallvec![0])) } global_phase += c2r.global_phase; let c2l = unitary_to_gate_sequence_inner( @@ -1058,7 +1059,7 @@ impl TwoQubitWeylDecomposition { ) .unwrap(); for gate in c2l.gates { - gate_sequence.push((gate.0, gate.1, smallvec![1])) + gate_sequence.push((gate.0.name().to_string(), gate.1, smallvec![1])) } global_phase += c2l.global_phase; self.weyl_gate( @@ -1077,7 +1078,7 @@ impl TwoQubitWeylDecomposition { ) .unwrap(); for gate in c1r.gates { - gate_sequence.push((gate.0, gate.1, smallvec![0])) + gate_sequence.push((gate.0.name().to_string(), gate.1, smallvec![0])) } global_phase += c2r.global_phase; let c1l = unitary_to_gate_sequence_inner( @@ -1090,7 +1091,7 @@ impl TwoQubitWeylDecomposition { ) .unwrap(); for gate in c1l.gates { - gate_sequence.push((gate.0, gate.1, smallvec![1])) + gate_sequence.push((gate.0.name().to_string(), gate.1, smallvec![1])) } Ok(TwoQubitGateSequence { gates: gate_sequence, @@ -1459,7 +1460,7 @@ impl TwoQubitBasisDecomposer { if let Some(sequence) = sequence { *global_phase += sequence.global_phase; for gate in sequence.gates { - gates.push((gate.0, gate.1, smallvec![qubit])); + gates.push((gate.0.name().to_string(), gate.1, smallvec![qubit])); } } } @@ -1847,13 +1848,13 @@ impl TwoQubitBasisDecomposer { for i in 0..best_nbasis as usize { if let Some(euler_decomp) = &euler_decompositions[2 * i] { for gate in &euler_decomp.gates { - gates.push((gate.0.clone(), gate.1.clone(), smallvec![0])); + gates.push((gate.0.name().to_string(), gate.1.clone(), smallvec![0])); } global_phase += euler_decomp.global_phase } if let Some(euler_decomp) = &euler_decompositions[2 * i + 1] { for gate in &euler_decomp.gates { - gates.push((gate.0.clone(), gate.1.clone(), smallvec![1])); + gates.push((gate.0.name().to_string(), gate.1.clone(), smallvec![1])); } global_phase += euler_decomp.global_phase } @@ -1861,13 +1862,13 @@ impl TwoQubitBasisDecomposer { } if let Some(euler_decomp) = &euler_decompositions[2 * best_nbasis as usize] { for gate in &euler_decomp.gates { - gates.push((gate.0.clone(), gate.1.clone(), smallvec![0])); + gates.push((gate.0.name().to_string(), gate.1.clone(), smallvec![0])); } global_phase += euler_decomp.global_phase } if let Some(euler_decomp) = &euler_decompositions[2 * best_nbasis as usize + 1] { for gate in &euler_decomp.gates { - gates.push((gate.0.clone(), gate.1.clone(), smallvec![1])); + gates.push((gate.0.name().to_string(), gate.1.clone(), smallvec![1])); } global_phase += euler_decomp.global_phase } diff --git a/crates/circuit/src/dag_node.rs b/crates/circuit/src/dag_node.rs index ffd7920a36fd..44ba5f7a6bf0 100644 --- a/crates/circuit/src/dag_node.rs +++ b/crates/circuit/src/dag_node.rs @@ -18,7 +18,7 @@ use crate::operations::Operation; use numpy::IntoPyArray; use pyo3::prelude::*; use pyo3::types::{PyDict, PyList, PySequence, PyString, PyTuple}; -use pyo3::{intern, IntoPy, PyObject, PyResult}; +use pyo3::{intern, IntoPy, PyObject, PyResult, ToPyObject}; use smallvec::smallvec; /// Parent class for DAGOpNode, DAGInNode, and DAGOutNode. @@ -135,6 +135,50 @@ impl DAGOpNode { )) } + #[staticmethod] + fn from_instruction( + py: Python, + instruction: CircuitInstruction, + dag: Option<&Bound>, + ) -> PyResult { + let qargs = instruction.qubits.clone_ref(py).into_bound(py); + let cargs = instruction.clbits.clone_ref(py).into_bound(py); + + let sort_key = match dag { + Some(dag) => { + let cache = dag + .getattr(intern!(py, "_key_cache"))? + .downcast_into_exact::()?; + let cache_key = PyTuple::new_bound(py, [&qargs, &cargs]); + match cache.get_item(&cache_key)? { + Some(key) => key, + None => { + let indices: PyResult> = qargs + .iter() + .chain(cargs.iter()) + .map(|bit| { + dag.call_method1(intern!(py, "find_bit"), (bit,))? + .getattr(intern!(py, "index")) + }) + .collect(); + let index_strs: Vec<_> = + indices?.into_iter().map(|i| format!("{:04}", i)).collect(); + let key = PyString::new_bound(py, index_strs.join(",").as_str()); + cache.set_item(&cache_key, &key)?; + key.into_any() + } + } + } + None => qargs.str()?.into_any(), + }; + let base = PyClassInitializer::from(DAGNode { _node_id: -1 }); + let sub = base.add_subclass(DAGOpNode { + instruction, + sort_key: sort_key.unbind(), + }); + Ok(Py::new(py, sub)?.to_object(py)) + } + fn __reduce__(slf: PyRef, py: Python) -> PyResult { let state = (slf.as_ref()._node_id, &slf.sort_key); Ok(( @@ -206,8 +250,8 @@ impl DAGOpNode { /// Returns the Instruction name corresponding to the op for this node #[getter] - fn get_name(&self, py: Python) -> PyObject { - self.instruction.operation.name().to_object(py) + fn get_name(&self) -> &str { + self.instruction.operation.name() } #[getter] diff --git a/crates/circuit/src/operations.rs b/crates/circuit/src/operations.rs index 8935e72e0ad5..77458e2fc35f 100644 --- a/crates/circuit/src/operations.rs +++ b/crates/circuit/src/operations.rs @@ -242,6 +242,12 @@ pub enum StandardGate { RZXGate = 52, } +impl ToPyObject for StandardGate { + fn to_object(&self, py: Python) -> PyObject { + self.into_py(py) + } +} + // TODO: replace all 34s (placeholders) with actual number static STANDARD_GATE_NUM_QUBITS: [u32; STANDARD_GATE_SIZE] = [ 1, 1, 1, 2, 2, 2, 3, 1, 1, 1, // 0-9 diff --git a/qiskit/dagcircuit/dagcircuit.py b/qiskit/dagcircuit/dagcircuit.py index a0d2c42be17e..b93a90e47f7b 100644 --- a/qiskit/dagcircuit/dagcircuit.py +++ b/qiskit/dagcircuit/dagcircuit.py @@ -54,6 +54,7 @@ from qiskit.dagcircuit.dagnode import DAGNode, DAGOpNode, DAGInNode, DAGOutNode from qiskit.circuit.bit import Bit from qiskit.pulse import Schedule +from qiskit._accelerate.euler_one_qubit_decomposer import collect_1q_runs_filter BitLocations = namedtuple("BitLocations", ("index", "registers")) # The allowable arguments to :meth:`DAGCircuit.copy_empty_like`'s ``vars_mode``. @@ -642,17 +643,17 @@ def _check_wires(self, args: Iterable[Bit | expr.Var], amap: dict[Bit | expr.Var if wire not in amap: raise DAGCircuitError(f"wire {wire} not found in {amap}") - def _increment_op(self, op): - if op.name in self._op_names: - self._op_names[op.name] += 1 + def _increment_op(self, op_name): + if op_name in self._op_names: + self._op_names[op_name] += 1 else: - self._op_names[op.name] = 1 + self._op_names[op_name] = 1 - def _decrement_op(self, op): - if self._op_names[op.name] == 1: - del self._op_names[op.name] + def _decrement_op(self, op_name): + if self._op_names[op_name] == 1: + del self._op_names[op_name] else: - self._op_names[op.name] -= 1 + self._op_names[op_name] -= 1 def copy_empty_like(self, *, vars_mode: _VarsMode = "alike"): """Return a copy of self with the same structure but empty. @@ -724,7 +725,7 @@ def _apply_op_node_back(self, node: DAGOpNode): additional = set(_additional_wires(node)).difference(node.cargs) node._node_id = self._multi_graph.add_node(node) - self._increment_op(node) + self._increment_op(node.name) # Add new in-edges from predecessors of the output nodes to the # operation node while deleting the old in-edges of the output nodes @@ -780,7 +781,7 @@ def apply_operation_back( node = DAGOpNode(op=op, qargs=qargs, cargs=cargs, dag=self) node._node_id = self._multi_graph.add_node(node) - self._increment_op(op) + self._increment_op(op.name) # Add new in-edges from predecessors of the output nodes to the # operation node while deleting the old in-edges of the output nodes @@ -832,7 +833,7 @@ def apply_operation_front( node = DAGOpNode(op=op, qargs=qargs, cargs=cargs, dag=self) node._node_id = self._multi_graph.add_node(node) - self._increment_op(op) + self._increment_op(op.name) # Add new out-edges to successors of the input nodes from the # operation node while deleting the old out-edges of the input nodes @@ -1379,10 +1380,10 @@ def replace_block_with_op( "Replacing the specified node block would introduce a cycle" ) from ex - self._increment_op(op) + self._increment_op(op.name) for nd in node_block: - self._decrement_op(nd.op) + self._decrement_op(nd.name) return new_node @@ -1593,7 +1594,7 @@ def edge_weight_map(wire): node_map = self._multi_graph.substitute_node_with_subgraph( node._node_id, in_dag._multi_graph, edge_map_fn, filter_fn, edge_weight_map ) - self._decrement_op(node.op) + self._decrement_op(node.name) variable_mapper = _classical_resource_map.VariableMapper( self.cregs.values(), wire_map, add_register=self.add_creg @@ -1624,7 +1625,7 @@ def edge_weight_map(wire): new_node = DAGOpNode(m_op, qargs=m_qargs, cargs=m_cargs, dag=self) new_node._node_id = new_node_index self._multi_graph[new_node_index] = new_node - self._increment_op(new_node.op) + self._increment_op(new_node.name) return {k: self._multi_graph[v] for k, v in node_map.items()} @@ -1696,17 +1697,17 @@ def substitute_node(self, node: DAGOpNode, op, inplace: bool = False, propagate_ if inplace: if op.name != node.op.name: - self._increment_op(op) - self._decrement_op(node.op) + self._increment_op(op.name) + self._decrement_op(node.name) node.op = op return node new_node = copy.copy(node) new_node.op = op self._multi_graph[node._node_id] = new_node - if op.name != node.op.name: - self._increment_op(op) - self._decrement_op(node.op) + if op.name != node.name: + self._increment_op(op.name) + self._decrement_op(node.name) return new_node def separable_circuits( @@ -1987,7 +1988,7 @@ def remove_op_node(self, node): self._multi_graph.remove_node_retain_edges( node._node_id, use_outgoing=False, condition=lambda edge1, edge2: edge1 == edge2 ) - self._decrement_op(node.op) + self._decrement_op(node.name) def remove_ancestors_of(self, node): """Remove all of the ancestor operation nodes of node.""" @@ -2152,19 +2153,7 @@ def filter_fn(node): def collect_1q_runs(self) -> list[list[DAGOpNode]]: """Return a set of non-conditional runs of 1q "op" nodes.""" - - def filter_fn(node): - return ( - isinstance(node, DAGOpNode) - and len(node.qargs) == 1 - and len(node.cargs) == 0 - and isinstance(node.op, Gate) - and hasattr(node.op, "__array__") - and getattr(node.op, "condition", None) is None - and not node.op.is_parameterized() - ) - - return rx.collect_runs(self._multi_graph, filter_fn) + return rx.collect_runs(self._multi_graph, collect_1q_runs_filter) def collect_2q_runs(self): """Return a set of non-conditional runs of 2q "op" nodes.""" diff --git a/qiskit/synthesis/one_qubit/one_qubit_decompose.py b/qiskit/synthesis/one_qubit/one_qubit_decompose.py index c84db761b7f0..f60f20f9524e 100644 --- a/qiskit/synthesis/one_qubit/one_qubit_decompose.py +++ b/qiskit/synthesis/one_qubit/one_qubit_decompose.py @@ -161,29 +161,19 @@ def build_circuit(self, gates, global_phase) -> QuantumCircuit | DAGCircuit: if len(gates) > 0 and isinstance(gates[0], tuple): lookup_gate = True - if self.use_dag: - from qiskit.dagcircuit import dagcircuit - - dag = dagcircuit.DAGCircuit() - dag.global_phase = global_phase - dag.add_qubits(qr) - for gate_entry in gates: - if lookup_gate: - gate = NAME_MAP[gate_entry[0]](*gate_entry[1]) - else: - gate = gate_entry - - dag.apply_operation_back(gate, (qr[0],), check=False) - return dag - else: - circuit = QuantumCircuit(qr, global_phase=global_phase) - for gate_entry in gates: - if lookup_gate: - gate = NAME_MAP[gate_entry[0]](*gate_entry[1]) - else: - gate = gate_entry - circuit._append(gate, [qr[0]], []) - return circuit + from qiskit.dagcircuit import dagcircuit + + dag = dagcircuit.DAGCircuit() + dag.global_phase = global_phase + dag.add_qubits(qr) + for gate_entry in gates: + if lookup_gate: + gate = NAME_MAP[gate_entry[0].name](*gate_entry[1]) + else: + gate = gate_entry.name + + dag.apply_operation_back(gate, (qr[0],), check=False) + return dag def __call__( self, @@ -225,11 +215,17 @@ def __call__( return self._decompose(unitary, simplify=simplify, atol=atol) def _decompose(self, unitary, simplify=True, atol=DEFAULT_ATOL): - circuit_sequence = euler_one_qubit_decomposer.unitary_to_gate_sequence( - unitary, [self.basis], 0, None, simplify, atol + if self.use_dag: + circuit_sequence = euler_one_qubit_decomposer.unitary_to_gate_sequence( + unitary, [self.basis], 0, None, simplify, atol + ) + circuit = self.build_circuit(circuit_sequence, circuit_sequence.global_phase) + return circuit + return QuantumCircuit._from_circuit_data( + euler_one_qubit_decomposer.unitary_to_circuit( + unitary, [self.basis], 0, None, simplify, atol + ) ) - circuit = self.build_circuit(circuit_sequence, circuit_sequence.global_phase) - return circuit @property def basis(self): diff --git a/qiskit/transpiler/passes/optimization/optimize_1q_decomposition.py b/qiskit/transpiler/passes/optimization/optimize_1q_decomposition.py index 3f8d07839c0f..04d95312aa6d 100644 --- a/qiskit/transpiler/passes/optimization/optimize_1q_decomposition.py +++ b/qiskit/transpiler/passes/optimization/optimize_1q_decomposition.py @@ -33,6 +33,7 @@ XGate, ) from qiskit.circuit import Qubit +from qiskit.circuit.quantumcircuitdata import CircuitInstruction from qiskit.dagcircuit.dagcircuit import DAGCircuit from qiskit.dagcircuit.dagnode import DAGOpNode @@ -110,16 +111,7 @@ def _build_error_map(self): else: return None - def _resynthesize_run(self, matrix, qubit=None): - """ - Re-synthesizes one 2x2 `matrix`, typically extracted via `dag.collect_1q_runs`. - - Returns the newly synthesized circuit in the indicated basis, or None - if no synthesis routine applied. - - When multiple synthesis options are available, it prefers the one with the lowest - error when the circuit is applied to `qubit`. - """ + def _get_decomposer(self, qubit=None): # include path for when target exists but target.num_qubits is None (BasicSimulator) if self._target is not None and self._target.num_qubits is not None: if qubit is not None: @@ -133,6 +125,19 @@ def _resynthesize_run(self, matrix, qubit=None): decomposers = _possible_decomposers(available_1q_basis) else: decomposers = self._global_decomposers + return decomposers + + def _resynthesize_run(self, matrix, qubit=None): + """ + Re-synthesizes one 2x2 `matrix`, typically extracted via `dag.collect_1q_runs`. + + Returns the newly synthesized circuit in the indicated basis, or None + if no synthesis routine applied. + + When multiple synthesis options are available, it prefers the one with the lowest + error when the circuit is applied to `qubit`. + """ + decomposers = self._get_decomposer(qubit) best_synth_circuit = euler_one_qubit_decomposer.unitary_to_gate_sequence( matrix, @@ -149,10 +154,13 @@ def _gate_sequence_to_dag(self, best_synth_circuit): out_dag.global_phase = best_synth_circuit.global_phase for gate_name, angles in best_synth_circuit: - out_dag.apply_operation_back(NAME_MAP[gate_name](*angles), qubits, check=False) + op = CircuitInstruction(gate_name, qubits=qubits, params=angles) + out_dag.apply_operation_back(op.operation, qubits, check=False) return out_dag - def _substitution_checks(self, dag, old_run, new_circ, basis, qubit): + def _substitution_checks( + self, dag, old_run, new_circ, basis, qubit, old_error=None, new_error=None + ): """ Returns `True` when it is recommended to replace `old_run` with `new_circ` over `basis`. """ @@ -176,11 +184,14 @@ def _substitution_checks(self, dag, old_run, new_circ, basis, qubit): # if we're outside of the basis set, we're obligated to logically decompose. # if we're outside of the set of gates for which we have physical definitions, # then we _try_ to decompose, using the results if we see improvement. - new_error = 0.0 - old_error = 0.0 if not uncalibrated_and_not_basis_p: - new_error = self._error(new_circ, qubit) - old_error = self._error(old_run, qubit) + if new_error is None: + new_error = self._error(new_circ, qubit) + if old_error is None: + old_error = self._error(old_run, qubit) + else: + new_error = 0.0 + old_error = 0.0 return ( uncalibrated_and_not_basis_p @@ -198,32 +209,47 @@ def run(self, dag): Returns: DAGCircuit: the optimized DAG. """ - runs = dag.collect_1q_runs() - for run in runs: + runs = [] + qubits = [] + bases = [] + for run in dag.collect_1q_runs(): qubit = dag.find_bit(run[0].qargs[0]).index - operator = run[0].op.to_matrix() - for node in run[1:]: - operator = node.op.to_matrix().dot(operator) - best_circuit_sequence = self._resynthesize_run(operator, qubit) - + 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 + ) + 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 and self._substitution_checks( - dag, run, best_circuit_sequence, basis, qubit - ): - for gate_name, angles in best_circuit_sequence: - op = NAME_MAP[gate_name](*angles) - node = DAGOpNode(NAME_MAP[gate_name](*angles), run[0].qargs, dag=dag) - node._node_id = dag._multi_graph.add_node(node) - dag._increment_op(op) - dag._multi_graph.insert_node_on_in_edges(node._node_id, run[0]._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) + 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(gate, qubits=qubit, params=angles) + node = DAGOpNode.from_instruction(op, dag=dag) + node._node_id = dag._multi_graph.add_node(node) + dag._increment_op(gate.name) + dag._multi_graph.insert_node_on_in_edges(node._node_id, 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 @@ -241,10 +267,7 @@ def _error(self, circuit, qubit): circuit, qubit, self.error_map ) else: - circuit_list = [(x.op.name, []) for x in circuit] - return euler_one_qubit_decomposer.compute_error_list( - circuit_list, qubit, self.error_map - ) + return euler_one_qubit_decomposer.compute_error_list(circuit, qubit, self.error_map) def _possible_decomposers(basis_set):