quantrs2_sim/
pauli.rs

1//! Pauli string evolution and operations.
2//!
3//! This module provides efficient operations for Pauli strings, including:
4//! - Pauli string construction and manipulation
5//! - Time evolution of Pauli observables
6//! - Commutation relations and algebra
7//! - Measurement expectation values
8
9use crate::prelude::SimulatorError;
10use scirs2_core::ndarray::Array2;
11use scirs2_core::parallel_ops::*;
12use scirs2_core::Complex64;
13use std::collections::HashMap;
14use std::fmt;
15
16use crate::error::Result;
17use crate::trotter::{DynamicCircuit, Hamiltonian, TrotterDecomposer, TrotterMethod};
18use quantrs2_core::gate::{multi::*, single::*, GateOp};
19use quantrs2_core::qubit::QubitId;
20
21/// Single Pauli operator type
22#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
23pub enum PauliOperator {
24    /// Identity operator
25    I,
26    /// Pauli-X operator
27    X,
28    /// Pauli-Y operator
29    Y,
30    /// Pauli-Z operator
31    Z,
32}
33
34impl fmt::Display for PauliOperator {
35    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
36        match self {
37            Self::I => write!(f, "I"),
38            Self::X => write!(f, "X"),
39            Self::Y => write!(f, "Y"),
40            Self::Z => write!(f, "Z"),
41        }
42    }
43}
44
45impl PauliOperator {
46    /// Parse from string
47    pub fn from_str(s: &str) -> Result<Self> {
48        match s.to_uppercase().as_str() {
49            "I" => Ok(Self::I),
50            "X" => Ok(Self::X),
51            "Y" => Ok(Self::Y),
52            "Z" => Ok(Self::Z),
53            _ => Err(SimulatorError::InvalidInput(format!(
54                "Invalid Pauli operator: {s}"
55            ))),
56        }
57    }
58
59    /// Get matrix representation
60    pub fn matrix(&self) -> Array2<Complex64> {
61        match self {
62            Self::I => Array2::from_shape_vec(
63                (2, 2),
64                vec![
65                    Complex64::new(1.0, 0.0),
66                    Complex64::new(0.0, 0.0),
67                    Complex64::new(0.0, 0.0),
68                    Complex64::new(1.0, 0.0),
69                ],
70            )
71            .unwrap(),
72            Self::X => Array2::from_shape_vec(
73                (2, 2),
74                vec![
75                    Complex64::new(0.0, 0.0),
76                    Complex64::new(1.0, 0.0),
77                    Complex64::new(1.0, 0.0),
78                    Complex64::new(0.0, 0.0),
79                ],
80            )
81            .unwrap(),
82            Self::Y => Array2::from_shape_vec(
83                (2, 2),
84                vec![
85                    Complex64::new(0.0, 0.0),
86                    Complex64::new(0.0, -1.0),
87                    Complex64::new(0.0, 1.0),
88                    Complex64::new(0.0, 0.0),
89                ],
90            )
91            .unwrap(),
92            Self::Z => Array2::from_shape_vec(
93                (2, 2),
94                vec![
95                    Complex64::new(1.0, 0.0),
96                    Complex64::new(0.0, 0.0),
97                    Complex64::new(0.0, 0.0),
98                    Complex64::new(-1.0, 0.0),
99                ],
100            )
101            .unwrap(),
102        }
103    }
104
105    /// Check if commutes with another Pauli
106    pub fn commutes_with(&self, other: &Self) -> bool {
107        match (self, other) {
108            (Self::I, _) | (_, Self::I) => true,
109            (a, b) if a == b => true,
110            _ => false,
111        }
112    }
113
114    /// Multiplication of Pauli operators (returns (result, phase))
115    pub const fn multiply(&self, other: &Self) -> (Self, Complex64) {
116        match (self, other) {
117            (Self::I, p) | (p, Self::I) => (*p, Complex64::new(1.0, 0.0)),
118            (Self::X, Self::X) | (Self::Y, Self::Y) | (Self::Z, Self::Z) => {
119                (Self::I, Complex64::new(1.0, 0.0))
120            }
121            (Self::X, Self::Y) => (Self::Z, Complex64::new(0.0, 1.0)),
122            (Self::Y, Self::X) => (Self::Z, Complex64::new(0.0, -1.0)),
123            (Self::Y, Self::Z) => (Self::X, Complex64::new(0.0, 1.0)),
124            (Self::Z, Self::Y) => (Self::X, Complex64::new(0.0, -1.0)),
125            (Self::Z, Self::X) => (Self::Y, Complex64::new(0.0, 1.0)),
126            (Self::X, Self::Z) => (Self::Y, Complex64::new(0.0, -1.0)),
127        }
128    }
129}
130
131/// A Pauli string is a tensor product of Pauli operators
132#[derive(Debug, Clone)]
133pub struct PauliString {
134    /// Pauli operators at each qubit position
135    pub operators: Vec<PauliOperator>,
136    /// Overall coefficient
137    pub coefficient: Complex64,
138    /// Number of qubits
139    pub num_qubits: usize,
140}
141
142impl PauliString {
143    /// Create a new Pauli string
144    pub fn new(num_qubits: usize) -> Self {
145        Self {
146            operators: vec![PauliOperator::I; num_qubits],
147            coefficient: Complex64::new(1.0, 0.0),
148            num_qubits,
149        }
150    }
151
152    /// Create from string representation like "XYZI"
153    pub fn from_string(pauli_str: &str, coefficient: Complex64) -> Result<Self> {
154        let operators: Result<Vec<_>> = pauli_str
155            .chars()
156            .map(|c| PauliOperator::from_str(&c.to_string()))
157            .collect();
158
159        Ok(Self {
160            operators: operators?,
161            coefficient,
162            num_qubits: pauli_str.len(),
163        })
164    }
165
166    /// Create from qubit indices and Pauli operators
167    pub fn from_ops(
168        num_qubits: usize,
169        ops: &[(usize, PauliOperator)],
170        coefficient: Complex64,
171    ) -> Result<Self> {
172        let mut pauli_string = Self::new(num_qubits);
173        pauli_string.coefficient = coefficient;
174
175        for &(qubit, op) in ops {
176            if qubit >= num_qubits {
177                return Err(SimulatorError::IndexOutOfBounds(qubit));
178            }
179            pauli_string.operators[qubit] = op;
180        }
181
182        Ok(pauli_string)
183    }
184
185    /// Set operator at specific qubit
186    pub fn set_operator(&mut self, qubit: usize, op: PauliOperator) -> Result<()> {
187        if qubit >= self.num_qubits {
188            return Err(SimulatorError::IndexOutOfBounds(qubit));
189        }
190        self.operators[qubit] = op;
191        Ok(())
192    }
193
194    /// Get operator at specific qubit
195    pub fn get_operator(&self, qubit: usize) -> Result<PauliOperator> {
196        if qubit >= self.num_qubits {
197            return Err(SimulatorError::IndexOutOfBounds(qubit));
198        }
199        Ok(self.operators[qubit])
200    }
201
202    /// Get non-identity operators
203    pub fn non_identity_ops(&self) -> Vec<(usize, PauliOperator)> {
204        self.operators
205            .iter()
206            .enumerate()
207            .filter(|(_, &op)| op != PauliOperator::I)
208            .map(|(i, &op)| (i, op))
209            .collect()
210    }
211
212    /// Check if this Pauli string commutes with another
213    pub fn commutes_with(&self, other: &Self) -> bool {
214        if self.num_qubits != other.num_qubits {
215            return false;
216        }
217
218        let mut anti_commute_count = 0;
219        for i in 0..self.num_qubits {
220            if !self.operators[i].commutes_with(&other.operators[i]) {
221                anti_commute_count += 1;
222            }
223        }
224
225        // Pauli strings commute if they anti-commute at an even number of positions
226        anti_commute_count % 2 == 0
227    }
228
229    /// Multiply two Pauli strings
230    pub fn multiply(&self, other: &Self) -> Result<Self> {
231        if self.num_qubits != other.num_qubits {
232            return Err(SimulatorError::DimensionMismatch(format!(
233                "Pauli strings have different lengths: {} vs {}",
234                self.num_qubits, other.num_qubits
235            )));
236        }
237
238        let mut result = Self::new(self.num_qubits);
239        let mut total_phase = self.coefficient * other.coefficient;
240
241        for i in 0..self.num_qubits {
242            let (op, phase) = self.operators[i].multiply(&other.operators[i]);
243            result.operators[i] = op;
244            total_phase *= phase;
245        }
246
247        result.coefficient = total_phase;
248        Ok(result)
249    }
250
251    /// Get weight (number of non-identity operators)
252    pub fn weight(&self) -> usize {
253        self.operators
254            .iter()
255            .filter(|&&op| op != PauliOperator::I)
256            .count()
257    }
258
259    /// Convert to Pauli string representation
260    pub fn pauli_string(&self) -> String {
261        self.operators.iter().map(|op| op.to_string()).collect()
262    }
263
264    /// Create time evolution circuit for this Pauli string
265    pub fn evolution_circuit(&self, time: f64) -> Result<DynamicCircuit> {
266        let mut circuit = DynamicCircuit::new(self.num_qubits);
267
268        if self.weight() == 0 {
269            // Identity operator - no gates needed
270            return Ok(circuit);
271        }
272
273        let non_identity = self.non_identity_ops();
274        let angle = -2.0 * self.coefficient.re * time;
275
276        if non_identity.len() == 1 {
277            // Single-qubit Pauli evolution
278            let (qubit, op) = non_identity[0];
279            match op {
280                PauliOperator::X => circuit.add_gate(Box::new(RotationX {
281                    target: QubitId::new(qubit as u32),
282                    theta: angle,
283                }))?,
284                PauliOperator::Y => circuit.add_gate(Box::new(RotationY {
285                    target: QubitId::new(qubit as u32),
286                    theta: angle,
287                }))?,
288                PauliOperator::Z => circuit.add_gate(Box::new(RotationZ {
289                    target: QubitId::new(qubit as u32),
290                    theta: angle,
291                }))?,
292                PauliOperator::I => {} // Should not happen
293            }
294        } else {
295            // Multi-qubit Pauli string evolution
296
297            // Apply basis rotations to convert all non-identity operators to Z
298            for &(qubit, op) in &non_identity {
299                match op {
300                    PauliOperator::X => circuit.add_gate(Box::new(Hadamard {
301                        target: QubitId::new(qubit as u32),
302                    }))?,
303                    PauliOperator::Y => {
304                        circuit.add_gate(Box::new(Hadamard {
305                            target: QubitId::new(qubit as u32),
306                        }))?;
307                        circuit.add_gate(Box::new(Phase {
308                            target: QubitId::new(qubit as u32),
309                        }))?;
310                    }
311                    PauliOperator::Z => {} // No basis change needed
312                    PauliOperator::I => {} // Should not happen
313                }
314            }
315
316            // Apply CNOT ladder to disentangle all Z operators to the last qubit
317            for i in 0..non_identity.len() - 1 {
318                circuit.add_gate(Box::new(CNOT {
319                    control: QubitId::new(non_identity[i].0 as u32),
320                    target: QubitId::new(non_identity[i + 1].0 as u32),
321                }))?;
322            }
323
324            // Apply Z rotation on the last qubit
325            circuit.add_gate(Box::new(RotationZ {
326                target: QubitId::new(non_identity[non_identity.len() - 1].0 as u32),
327                theta: angle,
328            }))?;
329
330            // Reverse CNOT ladder
331            for i in (0..non_identity.len() - 1).rev() {
332                circuit.add_gate(Box::new(CNOT {
333                    control: QubitId::new(non_identity[i].0 as u32),
334                    target: QubitId::new(non_identity[i + 1].0 as u32),
335                }))?;
336            }
337
338            // Reverse basis rotations
339            for &(qubit, op) in non_identity.iter().rev() {
340                match op {
341                    PauliOperator::X => circuit.add_gate(Box::new(Hadamard {
342                        target: QubitId::new(qubit as u32),
343                    }))?,
344                    PauliOperator::Y => {
345                        circuit.add_gate(Box::new(PhaseDagger {
346                            target: QubitId::new(qubit as u32),
347                        }))?;
348                        circuit.add_gate(Box::new(Hadamard {
349                            target: QubitId::new(qubit as u32),
350                        }))?;
351                    }
352                    PauliOperator::Z => {} // No basis change needed
353                    PauliOperator::I => {} // Should not happen
354                }
355            }
356        }
357
358        Ok(circuit)
359    }
360}
361
362impl fmt::Display for PauliString {
363    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
364        write!(f, "({}) {}", self.coefficient, self.pauli_string())
365    }
366}
367
368/// Sum of Pauli strings (Pauli operator)
369#[derive(Debug, Clone)]
370pub struct PauliOperatorSum {
371    /// Individual Pauli string terms
372    pub terms: Vec<PauliString>,
373    /// Number of qubits
374    pub num_qubits: usize,
375}
376
377impl PauliOperatorSum {
378    /// Create new empty sum
379    pub const fn new(num_qubits: usize) -> Self {
380        Self {
381            terms: Vec::new(),
382            num_qubits,
383        }
384    }
385
386    /// Add a Pauli string term
387    pub fn add_term(&mut self, pauli_string: PauliString) -> Result<()> {
388        if pauli_string.num_qubits != self.num_qubits {
389            return Err(SimulatorError::DimensionMismatch(format!(
390                "Pauli string has {} qubits, expected {}",
391                pauli_string.num_qubits, self.num_qubits
392            )));
393        }
394        self.terms.push(pauli_string);
395        Ok(())
396    }
397
398    /// Combine like terms
399    pub fn simplify(&mut self) {
400        let mut simplified_terms: HashMap<String, Complex64> = HashMap::new();
401
402        for term in &self.terms {
403            let key = format!("{term}");
404            *simplified_terms
405                .entry(key)
406                .or_insert(Complex64::new(0.0, 0.0)) += term.coefficient;
407        }
408
409        self.terms.clear();
410        for (pauli_str, coeff) in simplified_terms {
411            if coeff.norm() > 1e-15 {
412                if let Ok(term) = PauliString::from_string(&pauli_str, coeff) {
413                    self.terms.push(term);
414                }
415            }
416        }
417    }
418
419    /// Convert to Hamiltonian for Trotter evolution
420    pub fn to_hamiltonian(&self) -> Hamiltonian {
421        let mut ham = Hamiltonian::new(self.num_qubits);
422
423        for term in &self.terms {
424            let non_identity = term.non_identity_ops();
425
426            match non_identity.len() {
427                0 => {} // Identity term - ignore for Hamiltonian
428                1 => {
429                    let (qubit, op) = non_identity[0];
430                    let pauli_str = match op {
431                        PauliOperator::X => "X",
432                        PauliOperator::Y => "Y",
433                        PauliOperator::Z => "Z",
434                        PauliOperator::I => continue,
435                    };
436                    let _ = ham.add_single_pauli(qubit, pauli_str, term.coefficient.re);
437                }
438                2 => {
439                    let (q1, op1) = non_identity[0];
440                    let (q2, op2) = non_identity[1];
441                    let pauli1 = match op1 {
442                        PauliOperator::X => "X",
443                        PauliOperator::Y => "Y",
444                        PauliOperator::Z => "Z",
445                        PauliOperator::I => continue,
446                    };
447                    let pauli2 = match op2 {
448                        PauliOperator::X => "X",
449                        PauliOperator::Y => "Y",
450                        PauliOperator::Z => "Z",
451                        PauliOperator::I => continue,
452                    };
453                    let _ = ham.add_two_pauli(q1, q2, pauli1, pauli2, term.coefficient.re);
454                }
455                _ => {
456                    // Multi-qubit case
457                    let qubits: Vec<usize> = non_identity.iter().map(|&(q, _)| q).collect();
458                    let paulis: Vec<String> = non_identity
459                        .iter()
460                        .map(|&(_, op)| match op {
461                            PauliOperator::X => "X".to_string(),
462                            PauliOperator::Y => "Y".to_string(),
463                            PauliOperator::Z => "Z".to_string(),
464                            PauliOperator::I => "I".to_string(),
465                        })
466                        .collect();
467                    let _ = ham.add_pauli_string(qubits, paulis, term.coefficient.re);
468                }
469            }
470        }
471
472        ham
473    }
474
475    /// Time evolution using Trotter decomposition
476    pub fn time_evolution_circuit(
477        &self,
478        time: f64,
479        trotter_steps: usize,
480        method: TrotterMethod,
481    ) -> Result<DynamicCircuit> {
482        let hamiltonian = self.to_hamiltonian();
483        let decomposer = TrotterDecomposer::new(method, trotter_steps);
484        decomposer.decompose(&hamiltonian, time)
485    }
486
487    /// Direct time evolution without Trotter approximation (for single terms)
488    pub fn exact_evolution_circuit(&self, time: f64) -> Result<DynamicCircuit> {
489        let mut circuit = DynamicCircuit::new(self.num_qubits);
490
491        for term in &self.terms {
492            let term_circuit = term.evolution_circuit(time)?;
493            for gate in term_circuit.gates() {
494                circuit.add_gate(gate.clone())?;
495            }
496        }
497
498        Ok(circuit)
499    }
500}
501
502impl fmt::Display for PauliOperatorSum {
503    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
504        if self.terms.is_empty() {
505            write!(f, "0")
506        } else {
507            let term_strs: Vec<String> = self.terms.iter().map(|t| format!("{t}")).collect();
508            write!(f, "{}", term_strs.join(" + "))
509        }
510    }
511}
512
513/// Utilities for common Pauli string operations
514pub struct PauliUtils;
515
516impl PauliUtils {
517    /// Create a single-qubit Pauli string
518    pub fn single_qubit(
519        num_qubits: usize,
520        qubit: usize,
521        op: PauliOperator,
522        coeff: Complex64,
523    ) -> Result<PauliString> {
524        PauliString::from_ops(num_qubits, &[(qubit, op)], coeff)
525    }
526
527    /// Create an all-X Pauli string
528    pub fn all_x(num_qubits: usize, coeff: Complex64) -> PauliString {
529        let mut pauli = PauliString::new(num_qubits);
530        pauli.coefficient = coeff;
531        for i in 0..num_qubits {
532            pauli.operators[i] = PauliOperator::X;
533        }
534        pauli
535    }
536
537    /// Create an all-Z Pauli string
538    pub fn all_z(num_qubits: usize, coeff: Complex64) -> PauliString {
539        let mut pauli = PauliString::new(num_qubits);
540        pauli.coefficient = coeff;
541        for i in 0..num_qubits {
542            pauli.operators[i] = PauliOperator::Z;
543        }
544        pauli
545    }
546
547    /// Create random Pauli string
548    pub fn random(num_qubits: usize, weight: usize, coeff: Complex64) -> Result<PauliString> {
549        if weight > num_qubits {
550            return Err(SimulatorError::InvalidInput(
551                "Weight cannot exceed number of qubits".to_string(),
552            ));
553        }
554
555        let mut pauli = PauliString::new(num_qubits);
556        pauli.coefficient = coeff;
557
558        // Randomly select positions for non-identity operators
559        let mut positions: Vec<usize> = (0..num_qubits).collect();
560        fastrand::shuffle(&mut positions);
561
562        let ops = [PauliOperator::X, PauliOperator::Y, PauliOperator::Z];
563
564        for i in 0..weight {
565            pauli.operators[positions[i]] = ops[fastrand::usize(0..3)];
566        }
567
568        Ok(pauli)
569    }
570
571    /// Check if a set of Pauli strings are mutually commuting
572    pub fn are_mutually_commuting(pauli_strings: &[PauliString]) -> bool {
573        for i in 0..pauli_strings.len() {
574            for j in i + 1..pauli_strings.len() {
575                if !pauli_strings[i].commutes_with(&pauli_strings[j]) {
576                    return false;
577                }
578            }
579        }
580        true
581    }
582
583    /// Find maximal set of mutually commuting Pauli strings
584    pub fn maximal_commuting_set(pauli_strings: &[PauliString]) -> Vec<usize> {
585        let mut commuting_set = Vec::new();
586
587        for (i, pauli) in pauli_strings.iter().enumerate() {
588            let mut can_add = true;
589            for &j in &commuting_set {
590                if !pauli.commutes_with(&pauli_strings[j]) {
591                    can_add = false;
592                    break;
593                }
594            }
595
596            if can_add {
597                commuting_set.push(i);
598            }
599        }
600
601        commuting_set
602    }
603}
604
605#[cfg(test)]
606mod tests {
607    use super::*;
608
609    #[test]
610    fn test_pauli_operator_multiply() {
611        let (result, phase) = PauliOperator::X.multiply(&PauliOperator::Y);
612        assert_eq!(result, PauliOperator::Z);
613        assert_eq!(phase, Complex64::new(0.0, 1.0));
614
615        let (result, phase) = PauliOperator::Y.multiply(&PauliOperator::X);
616        assert_eq!(result, PauliOperator::Z);
617        assert_eq!(phase, Complex64::new(0.0, -1.0));
618    }
619
620    #[test]
621    fn test_pauli_string_creation() {
622        let pauli = PauliString::from_string("XYZ", Complex64::new(1.0, 0.0)).unwrap();
623        assert_eq!(pauli.num_qubits, 3);
624        assert_eq!(pauli.operators[0], PauliOperator::X);
625        assert_eq!(pauli.operators[1], PauliOperator::Y);
626        assert_eq!(pauli.operators[2], PauliOperator::Z);
627    }
628
629    #[test]
630    fn test_pauli_string_multiply() {
631        let p1 = PauliString::from_string("XY", Complex64::new(1.0, 0.0)).unwrap();
632        let p2 = PauliString::from_string("YZ", Complex64::new(1.0, 0.0)).unwrap();
633
634        let result = p1.multiply(&p2).unwrap();
635        assert_eq!(result.operators[0], PauliOperator::Z);
636        assert_eq!(result.operators[1], PauliOperator::X);
637        assert_eq!(result.coefficient, Complex64::new(-1.0, 0.0));
638    }
639
640    #[test]
641    fn test_pauli_string_commutation() {
642        let p1 = PauliString::from_string("XY", Complex64::new(1.0, 0.0)).unwrap();
643        let p2 = PauliString::from_string("ZI", Complex64::new(1.0, 0.0)).unwrap();
644        let p3 = PauliString::from_string("XI", Complex64::new(1.0, 0.0)).unwrap();
645
646        assert!(!p1.commutes_with(&p2)); // XY and ZI anti-commute at first qubit
647        assert!(p1.commutes_with(&p3)); // XY and XI commute
648    }
649
650    #[test]
651    fn test_pauli_string_weight() {
652        let p1 = PauliString::from_string("XIYZ", Complex64::new(1.0, 0.0)).unwrap();
653        assert_eq!(p1.weight(), 3);
654
655        let p2 = PauliString::from_string("IIII", Complex64::new(1.0, 0.0)).unwrap();
656        assert_eq!(p2.weight(), 0);
657    }
658
659    #[test]
660    fn test_pauli_operator_sum() {
661        let mut sum = PauliOperatorSum::new(2);
662
663        let p1 = PauliString::from_string("XX", Complex64::new(1.0, 0.0)).unwrap();
664        let p2 = PauliString::from_string("YY", Complex64::new(0.5, 0.0)).unwrap();
665
666        sum.add_term(p1).unwrap();
667        sum.add_term(p2).unwrap();
668
669        assert_eq!(sum.terms.len(), 2);
670    }
671
672    #[test]
673    fn test_evolution_circuit_single_qubit() {
674        let pauli = PauliString::from_string("X", Complex64::new(1.0, 0.0)).unwrap();
675        let circuit = pauli.evolution_circuit(1.0).unwrap();
676        assert!(circuit.gate_count() > 0);
677    }
678
679    #[test]
680    fn test_evolution_circuit_multi_qubit() {
681        let pauli = PauliString::from_string("XYZ", Complex64::new(1.0, 0.0)).unwrap();
682        let circuit = pauli.evolution_circuit(1.0).unwrap();
683        assert!(circuit.gate_count() > 0);
684    }
685
686    #[test]
687    fn test_utils_single_qubit() {
688        let pauli =
689            PauliUtils::single_qubit(3, 1, PauliOperator::X, Complex64::new(1.0, 0.0)).unwrap();
690        assert_eq!(pauli.operators[0], PauliOperator::I);
691        assert_eq!(pauli.operators[1], PauliOperator::X);
692        assert_eq!(pauli.operators[2], PauliOperator::I);
693    }
694
695    #[test]
696    fn test_mutually_commuting() {
697        let p1 = PauliString::from_string("XX", Complex64::new(1.0, 0.0)).unwrap();
698        let p2 = PauliString::from_string("ZZ", Complex64::new(1.0, 0.0)).unwrap();
699        let p3 = PauliString::from_string("YY", Complex64::new(1.0, 0.0)).unwrap();
700
701        let paulis = vec![p1, p2, p3];
702        assert!(PauliUtils::are_mutually_commuting(&paulis));
703    }
704}