[add]: add support for non-clifford gates

3 files changed, 452 insertions, 25 deletions

diff --git a/libpsi-core/src/core/circuit.rs b/libpsi-core/src/core/circuit.rs
index 5bb17ff..e5d7847 100644
--- a/libpsi-core/src/core/circuit.rs
+++ b/libpsi-core/src/core/circuit.rs
@@ -1,7 +1,7 @@
 use super::{CustomGate, QuantumState, Runtime};
 use crate::{format_amplitude, format_probability, Vector};
-use std::sync::Arc;
 use core::fmt;
+use std::sync::Arc;
 
 #[derive(Clone)]
 pub enum GateOp {
@@ -11,9 +11,24 @@ pub enum GateOp {
     Z(usize),
     S(usize),
     T(usize),
+    Sdg(usize),
+    Tdg(usize),
+    Sx(usize),
+    Sxdg(usize),
+    Rx(usize, f64),
+    Ry(usize, f64),
+    Rz(usize, f64),
+    P(usize, f64),
+    U1(usize, f64),
+    U2(usize, f64, f64),
+    U3(usize, f64, f64, f64),
     CNOT(usize, usize),
     CZ(usize, usize),
     SWAP(usize, usize),
+    CRx(usize, usize, f64),
+    CRy(usize, usize, f64),
+    CRz(usize, usize, f64),
+    CP(usize, usize, f64),
     CCNOT(usize, usize, usize),
     CSWAP(usize, usize, usize),
     Measure(usize, usize),
@@ -29,6 +44,21 @@ impl GateOp {
             GateOp::Z(_) => "Z",
             GateOp::S(_) => "S",
             GateOp::T(_) => "T",
+            GateOp::Sdg(_) => "S†",
+            GateOp::Tdg(_) => "T†",
+            GateOp::Sx(_) => "√X",
+            GateOp::Sxdg(_) => "√X†",
+            GateOp::Rx(_, _) => "Rx",
+            GateOp::Ry(_, _) => "Ry",
+            GateOp::Rz(_, _) => "Rz",
+            GateOp::P(_, _) => "P",
+            GateOp::U1(_, _) => "U1",
+            GateOp::U2(_, _, _) => "U2",
+            GateOp::U3(_, _, _, _) => "U3",
+            GateOp::CRx(_, _, _) => "CRx",
+            GateOp::CRy(_, _, _) => "CRy",
+            GateOp::CRz(_, _, _) => "CRz",
+            GateOp::CP(_, _, _) => "CP",
             GateOp::CNOT(_, _) => "CNOT",
             GateOp::CZ(_, _) => "CZ",
             GateOp::SWAP(_, _) => "SWAP",
@@ -41,8 +71,30 @@ impl GateOp {
 
     pub fn quantum_targets(&self) -> Vec<usize> {
         match self {
-            GateOp::H(t) | GateOp::X(t) | GateOp::Y(t) | GateOp::Z(t) | GateOp::S(t) | GateOp::T(t) => vec![*t],
-            GateOp::CNOT(c, t) | GateOp::CZ(c, t) | GateOp::SWAP(c, t) => vec![*c, *t],
+            GateOp::H(t)
+            | GateOp::X(t)
+            | GateOp::Y(t)
+            | GateOp::Z(t)
+            | GateOp::S(t)
+            | GateOp::T(t)
+            | GateOp::Sdg(t)
+            | GateOp::Tdg(t)
+            | GateOp::Sx(t)
+            | GateOp::Sxdg(t)
+            | GateOp::Rx(t, _)
+            | GateOp::Ry(t, _)
+            | GateOp::Rz(t, _)
+            | GateOp::P(t, _)
+            | GateOp::U1(t, _)
+            | GateOp::U2(t, _, _)
+            | GateOp::U3(t, _, _, _) => vec![*t],
+            GateOp::CNOT(c, t)
+            | GateOp::CZ(c, t)
+            | GateOp::SWAP(c, t)
+            | GateOp::CRx(c, t, _)
+            | GateOp::CRy(c, t, _)
+            | GateOp::CRz(c, t, _)
+            | GateOp::CP(c, t, _) => vec![*c, *t],
             GateOp::CCNOT(c1, c2, t) | GateOp::CSWAP(c1, c2, t) => vec![*c1, *c2, *t],
             GateOp::Measure(q, _) => vec![*q],
             GateOp::Custom(_, targets) => targets.clone(),
@@ -63,6 +115,27 @@ impl GateOp {
     pub fn is_custom(&self) -> bool {
         matches!(self, GateOp::Custom(_, _))
     }
+
+    pub fn is_non_clifford(&self) -> bool {
+        matches!(
+            self,
+            GateOp::T(_)
+                | GateOp::Tdg(_)
+                | GateOp::Sx(_)
+                | GateOp::Sxdg(_)
+                | GateOp::Rx(_, _)
+                | GateOp::Ry(_, _)
+                | GateOp::Rz(_, _)
+                | GateOp::P(_, _)
+                | GateOp::U1(_, _)
+                | GateOp::U2(_, _, _)
+                | GateOp::U3(_, _, _, _)
+                | GateOp::CRx(_, _, _)
+                | GateOp::CRy(_, _, _)
+                | GateOp::CRz(_, _, _)
+                | GateOp::CP(_, _, _)
+        )
+    }
 }
 
 pub struct QuantumCircuit {
@@ -164,6 +237,96 @@ impl QuantumCircuit {
         self
     }
 
+    pub fn sdg(&mut self, target: usize) -> &mut Self {
+        self.operations.push(GateOp::Sdg(target));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn tdg(&mut self, target: usize) -> &mut Self {
+        self.operations.push(GateOp::Tdg(target));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn sx(&mut self, target: usize) -> &mut Self {
+        self.operations.push(GateOp::Sx(target));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn sxdg(&mut self, target: usize) -> &mut Self {
+        self.operations.push(GateOp::Sxdg(target));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn rx(&mut self, target: usize, theta: f64) -> &mut Self {
+        self.operations.push(GateOp::Rx(target, theta));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn ry(&mut self, target: usize, theta: f64) -> &mut Self {
+        self.operations.push(GateOp::Ry(target, theta));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn rz(&mut self, target: usize, theta: f64) -> &mut Self {
+        self.operations.push(GateOp::Rz(target, theta));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn p(&mut self, target: usize, theta: f64) -> &mut Self {
+        self.operations.push(GateOp::P(target, theta));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn u1(&mut self, target: usize, lambda: f64) -> &mut Self {
+        self.operations.push(GateOp::U1(target, lambda));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn u2(&mut self, target: usize, phi: f64, lambda: f64) -> &mut Self {
+        self.operations.push(GateOp::U2(target, phi, lambda));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn u3(&mut self, target: usize, theta: f64, phi: f64, lambda: f64) -> &mut Self {
+        self.operations.push(GateOp::U3(target, theta, phi, lambda));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn crx(&mut self, control: usize, target: usize, theta: f64) -> &mut Self {
+        self.operations.push(GateOp::CRx(control, target, theta));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn cry(&mut self, control: usize, target: usize, theta: f64) -> &mut Self {
+        self.operations.push(GateOp::CRy(control, target, theta));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn crz(&mut self, control: usize, target: usize, theta: f64) -> &mut Self {
+        self.operations.push(GateOp::CRz(control, target, theta));
+        self.computed_state = None;
+        self
+    }
+
+    pub fn cp(&mut self, control: usize, target: usize, theta: f64) -> &mut Self {
+        self.operations.push(GateOp::CP(control, target, theta));
+        self.computed_state = None;
+        self
+    }
+
     pub fn cnot(&mut self, control: usize, target: usize) -> &mut Self {
         self.operations.push(GateOp::CNOT(control, target));
         self.computed_state = None;
@@ -187,7 +350,8 @@ impl QuantumCircuit {
     }
 
     pub fn ccnot(&mut self, control1: usize, control2: usize, target: usize) -> &mut Self {
-        self.operations.push(GateOp::CCNOT(control1, control2, target));
+        self.operations
+            .push(GateOp::CCNOT(control1, control2, target));
         self.computed_state = None;
         self
     }
@@ -197,7 +361,8 @@ impl QuantumCircuit {
     }
 
     pub fn cswap(&mut self, control: usize, target1: usize, target2: usize) -> &mut Self {
-        self.operations.push(GateOp::CSWAP(control, target1, target2));
+        self.operations
+            .push(GateOp::CSWAP(control, target1, target2));
         self.computed_state = None;
         self
     }
@@ -222,13 +387,15 @@ impl QuantumCircuit {
     }
 
     pub fn custom(&mut self, gate: &Arc<CustomGate>, targets: &[usize]) -> &mut Self {
-        self.operations.push(GateOp::Custom(Arc::clone(gate), targets.to_vec()));
+        self.operations
+            .push(GateOp::Custom(Arc::clone(gate), targets.to_vec()));
         self.computed_state = None;
         self
     }
 
     pub fn apply_custom(&mut self, gate: CustomGate, targets: &[usize]) -> &mut Self {
-        self.operations.push(GateOp::Custom(Arc::new(gate), targets.to_vec()));
+        self.operations
+            .push(GateOp::Custom(Arc::new(gate), targets.to_vec()));
         self.computed_state = None;
         self
     }
@@ -268,12 +435,18 @@ impl QuantumCircuit {
 
 impl fmt::Display for QuantumCircuit {
     fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
-        writeln!(f, "QuantumCircuit ({} qubits, {} classical)", self.num_qubits, self.num_classical)?;
+        writeln!(
+            f,
+            "QuantumCircuit ({} qubits, {} classical)",
+            self.num_qubits, self.num_classical
+        )?;
         writeln!(f, "Operations:")?;
         for (i, op) in self.operations.iter().enumerate() {
             match op {
                 GateOp::Measure(q, c) => writeln!(f, "  {}: {} q{} → c{}", i, op.name(), q, c)?,
-                GateOp::Custom(gate, targets) => writeln!(f, "  {}: [{}] on {:?}", i, gate.name, targets)?,
+                GateOp::Custom(gate, targets) => {
+                    writeln!(f, "  {}: [{}] on {:?}", i, gate.name, targets)?
+                }
                 _ => writeln!(f, "  {}: {} on {:?}", i, op.name(), op.quantum_targets())?,
             }
         }
diff --git a/libpsi-core/src/core/gates.rs b/libpsi-core/src/core/gates.rs
index 4f922a8..cee4ebc 100644
--- a/libpsi-core/src/core/gates.rs
+++ b/libpsi-core/src/core/gates.rs
@@ -1,4 +1,100 @@
-use crate::{complex, matrix, QuantumGate};
+use crate::{complex, matrix, Complex, Matrix, QuantumGate};
+use std::f64::consts::FRAC_1_SQRT_2;
+
+pub fn rx_matrix(theta: f64) -> Matrix<Complex<f64>> {
+    let cos = (theta / 2.0).cos();
+    let sin = (theta / 2.0).sin();
+    matrix!(
+        [complex!(cos, 0.0), complex!(0.0, -sin)];
+        [complex!(0.0, -sin), complex!(cos, 0.0)]
+    )
+}
+
+pub fn ry_matrix(theta: f64) -> Matrix<Complex<f64>> {
+    let cos = (theta / 2.0).cos();
+    let sin = (theta / 2.0).sin();
+    matrix!(
+        [complex!(cos, 0.0), complex!(-sin, 0.0)];
+        [complex!(sin, 0.0), complex!(cos, 0.0)]
+    )
+}
+
+pub fn rz_matrix(theta: f64) -> Matrix<Complex<f64>> {
+    let half = theta / 2.0;
+    matrix!(
+        [complex!(half.cos(), -half.sin()), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(half.cos(), half.sin())]
+    )
+}
+
+pub fn p_matrix(theta: f64) -> Matrix<Complex<f64>> {
+    matrix!(
+        [complex!(1.0, 0.0), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(theta.cos(), theta.sin())]
+    )
+}
+
+pub fn u1_matrix(lambda: f64) -> Matrix<Complex<f64>> {
+    p_matrix(lambda)
+}
+
+pub fn u2_matrix(phi: f64, lambda: f64) -> Matrix<Complex<f64>> {
+    let inv_sqrt2 = FRAC_1_SQRT_2;
+    matrix!(
+        [complex!(inv_sqrt2, 0.0), complex!(-inv_sqrt2 * lambda.cos(), -inv_sqrt2 * lambda.sin())];
+        [complex!(inv_sqrt2 * phi.cos(), inv_sqrt2 * phi.sin()), complex!((phi + lambda).cos() * inv_sqrt2, (phi + lambda).sin() * inv_sqrt2)]
+    )
+}
+
+pub fn u3_matrix(theta: f64, phi: f64, lambda: f64) -> Matrix<Complex<f64>> {
+    let cos = (theta / 2.0).cos();
+    let sin = (theta / 2.0).sin();
+    matrix!(
+        [complex!(cos, 0.0), complex!(-sin * lambda.cos(), -sin * lambda.sin())];
+        [complex!(sin * phi.cos(), sin * phi.sin()), complex!(cos * (phi + lambda).cos(), cos * (phi + lambda).sin())]
+    )
+}
+
+pub fn crx_matrix(theta: f64) -> Matrix<Complex<f64>> {
+    let cos = (theta / 2.0).cos();
+    let sin = (theta / 2.0).sin();
+    matrix!(
+        [complex!(1.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(1.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(cos, 0.0), complex!(0.0, -sin)];
+        [complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(0.0, -sin), complex!(cos, 0.0)]
+    )
+}
+
+pub fn cry_matrix(theta: f64) -> Matrix<Complex<f64>> {
+    let cos = (theta / 2.0).cos();
+    let sin = (theta / 2.0).sin();
+    matrix!(
+        [complex!(1.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(1.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(cos, 0.0), complex!(-sin, 0.0)];
+        [complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(sin, 0.0), complex!(cos, 0.0)]
+    )
+}
+
+pub fn crz_matrix(theta: f64) -> Matrix<Complex<f64>> {
+    let half = theta / 2.0;
+    matrix!(
+        [complex!(1.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(1.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(half.cos(), -half.sin()), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(half.cos(), half.sin())]
+    )
+}
+
+pub fn cp_matrix(theta: f64) -> Matrix<Complex<f64>> {
+    matrix!(
+        [complex!(1.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(1.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(1.0, 0.0), complex!(0.0, 0.0)];
+        [complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(0.0, 0.0), complex!(theta.cos(), theta.sin())]
+    )
+}
 
 #[rustfmt::skip]
 lazy_static::lazy_static! {
@@ -45,6 +141,34 @@ lazy_static::lazy_static! {
         num_qubits: 1,
     };
     
+    pub static ref SDG_GATE: QuantumGate<'static> = QuantumGate {
+        name: "S†",
+        matrix: matrix!([complex!(1.0, 0.0), complex!(0.0, 0.0)];
+                        [complex!(0.0, 0.0), complex!(0.0, -1.0)]),
+        num_qubits: 1,
+    };
+    
+    pub static ref TDG_GATE: QuantumGate<'static> = QuantumGate {
+        name: "T†",
+        matrix: matrix!([complex!(1.0, 0.0), complex!(0.0, 0.0)];
+                        [complex!(0.0, 0.0), complex!(core::f64::consts::FRAC_1_SQRT_2, -core::f64::consts::FRAC_1_SQRT_2)]),
+        num_qubits: 1,
+    };
+    
+    pub static ref SX_GATE: QuantumGate<'static> = QuantumGate {
+        name: "√X",
+        matrix: matrix!([complex!(0.5, 0.5), complex!(0.5, -0.5)];
+                        [complex!(0.5, -0.5), complex!(0.5, 0.5)]),
+        num_qubits: 1,
+    };
+    
+    pub static ref SXDG_GATE: QuantumGate<'static> = QuantumGate {
+        name: "√X†",
+        matrix: matrix!([complex!(0.5, -0.5), complex!(0.5, 0.5)];
+                        [complex!(0.5, 0.5), complex!(0.5, -0.5)]),
+        num_qubits: 1,
+    };
+    
     pub static ref IDENTITY: QuantumGate<'static> = QuantumGate {
         name: "I",
         matrix: matrix!([complex!(1.0, 0.0), complex!(0.0, 0.0)];
diff --git a/libpsi-core/src/core/runtime.rs b/libpsi-core/src/core/runtime.rs
index 7b73c47..0a78f0f 100644
--- a/libpsi-core/src/core/runtime.rs
+++ b/libpsi-core/src/core/runtime.rs
@@ -1,5 +1,10 @@
 use super::{GateOp, QuantumGate, QuantumRegister, QuantumState};
-use crate::gates::*;
+use crate::gates::{
+    cp_matrix, crx_matrix, cry_matrix, crz_matrix, p_matrix, rx_matrix, ry_matrix,
+    rz_matrix, u1_matrix, u2_matrix, u3_matrix, CNOT, CZ, FREDKIN, HADAMARD,
+    PAULI_X, PAULI_Y, PAULI_Z, SDG_GATE, SWAP, SXDG_GATE, SX_GATE, S_GATE, TDG_GATE,
+    TOFFOLI, T_GATE,
+};
 use crate::maths::vector::Vector;
 use crate::{complex, Complex, Matrix};
 use rayon::prelude::*;
@@ -48,17 +53,118 @@ impl Runtime {
 
         for op in operations {
             match op {
+                // Clifford gates
                 GateOp::H(t) => register.apply_gate(&HADAMARD, &[*t]),
                 GateOp::X(t) => register.apply_gate(&PAULI_X, &[*t]),
                 GateOp::Y(t) => register.apply_gate(&PAULI_Y, &[*t]),
                 GateOp::Z(t) => register.apply_gate(&PAULI_Z, &[*t]),
                 GateOp::S(t) => register.apply_gate(&S_GATE, &[*t]),
-                GateOp::T(t) => register.apply_gate(&T_GATE, &[*t]),
                 GateOp::CNOT(c, t) => register.apply_gate(&CNOT, &[*c, *t]),
                 GateOp::CZ(c, t) => register.apply_gate(&CZ, &[*c, *t]),
                 GateOp::SWAP(a, b) => register.apply_gate(&SWAP, &[*a, *b]),
                 GateOp::CCNOT(c1, c2, t) => register.apply_gate(&TOFFOLI, &[*c1, *c2, *t]),
                 GateOp::CSWAP(c, t1, t2) => register.apply_gate(&FREDKIN, &[*c, *t1, *t2]),
+                
+                // Non-Clifford fixed gates
+                GateOp::T(t) => register.apply_gate(&T_GATE, &[*t]),
+                GateOp::Sdg(t) => register.apply_gate(&SDG_GATE, &[*t]),
+                GateOp::Tdg(t) => register.apply_gate(&TDG_GATE, &[*t]),
+                GateOp::Sx(t) => register.apply_gate(&SX_GATE, &[*t]),
+                GateOp::Sxdg(t) => register.apply_gate(&SXDG_GATE, &[*t]),
+                
+                // Parametric single-qubit gates (non-Clifford for most angles)
+                GateOp::Rx(t, theta) => {
+                    let gate = QuantumGate {
+                        name: "Rx",
+                        matrix: rx_matrix(*theta),
+                        num_qubits: 1,
+                    };
+                    register.apply_gate(&gate, &[*t]);
+                }
+                GateOp::Ry(t, theta) => {
+                    let gate = QuantumGate {
+                        name: "Ry",
+                        matrix: ry_matrix(*theta),
+                        num_qubits: 1,
+                    };
+                    register.apply_gate(&gate, &[*t]);
+                }
+                GateOp::Rz(t, theta) => {
+                    let gate = QuantumGate {
+                        name: "Rz",
+                        matrix: rz_matrix(*theta),
+                        num_qubits: 1,
+                    };
+                    register.apply_gate(&gate, &[*t]);
+                }
+                GateOp::P(t, theta) => {
+                    let gate = QuantumGate {
+                        name: "P",
+                        matrix: p_matrix(*theta),
+                        num_qubits: 1,
+                    };
+                    register.apply_gate(&gate, &[*t]);
+                }
+                GateOp::U1(t, lambda) => {
+                    let gate = QuantumGate {
+                        name: "U1",
+                        matrix: u1_matrix(*lambda),
+                        num_qubits: 1,
+                    };
+                    register.apply_gate(&gate, &[*t]);
+                }
+                GateOp::U2(t, phi, lambda) => {
+                    let gate = QuantumGate {
+                        name: "U2",
+                        matrix: u2_matrix(*phi, *lambda),
+                        num_qubits: 1,
+                    };
+                    register.apply_gate(&gate, &[*t]);
+                }
+                GateOp::U3(t, theta, phi, lambda) => {
+                    let gate = QuantumGate {
+                        name: "U3",
+                        matrix: u3_matrix(*theta, *phi, *lambda),
+                        num_qubits: 1,
+                    };
+                    register.apply_gate(&gate, &[*t]);
+                }
+                
+                // Controlled parametric gates
+                GateOp::CRx(c, t, theta) => {
+                    let gate = QuantumGate {
+                        name: "CRx",
+                        matrix: crx_matrix(*theta),
+                        num_qubits: 2,
+                    };
+                    register.apply_gate(&gate, &[*c, *t]);
+                }
+                GateOp::CRy(c, t, theta) => {
+                    let gate = QuantumGate {
+                        name: "CRy",
+                        matrix: cry_matrix(*theta),
+                        num_qubits: 2,
+                    };
+                    register.apply_gate(&gate, &[*c, *t]);
+                }
+                GateOp::CRz(c, t, theta) => {
+                    let gate = QuantumGate {
+                        name: "CRz",
+                        matrix: crz_matrix(*theta),
+                        num_qubits: 2,
+                    };
+                    register.apply_gate(&gate, &[*c, *t]);
+                }
+                GateOp::CP(c, t, theta) => {
+                    let gate = QuantumGate {
+                        name: "CP",
+                        matrix: cp_matrix(*theta),
+                        num_qubits: 2,
+                    };
+                    register.apply_gate(&gate, &[*c, *t]);
+                }
+                
+                // Measurement and custom gates
                 GateOp::Measure(_, _) => {}
                 GateOp::Custom(gate, targets) => {
                     let quantum_gate = gate.to_quantum_gate();
@@ -83,18 +189,42 @@ impl Runtime {
         state[0] = complex!(1.0, 0.0);
 
         for op in operations {
-            let (gate, targets): (&QuantumGate, Vec<usize>) = match op {
-                GateOp::H(t) => (&HADAMARD, vec![*t]),
-                GateOp::X(t) => (&PAULI_X, vec![*t]),
-                GateOp::Y(t) => (&PAULI_Y, vec![*t]),
-                GateOp::Z(t) => (&PAULI_Z, vec![*t]),
-                GateOp::S(t) => (&S_GATE, vec![*t]),
-                GateOp::T(t) => (&T_GATE, vec![*t]),
-                GateOp::CNOT(c, t) => (&CNOT, vec![*c, *t]),
-                GateOp::CZ(c, t) => (&CZ, vec![*c, *t]),
-                GateOp::SWAP(a, b) => (&SWAP, vec![*a, *b]),
-                GateOp::CCNOT(c1, c2, t) => (&TOFFOLI, vec![*c1, *c2, *t]),
-                GateOp::CSWAP(c, t1, t2) => (&FREDKIN, vec![*c, *t1, *t2]),
+            let (gate_matrix, targets): (Matrix<Complex<f64>>, Vec<usize>) = match op {
+                // Clifford gates
+                GateOp::H(t) => (HADAMARD.matrix.clone(), vec![*t]),
+                GateOp::X(t) => (PAULI_X.matrix.clone(), vec![*t]),
+                GateOp::Y(t) => (PAULI_Y.matrix.clone(), vec![*t]),
+                GateOp::Z(t) => (PAULI_Z.matrix.clone(), vec![*t]),
+                GateOp::S(t) => (S_GATE.matrix.clone(), vec![*t]),
+                GateOp::CNOT(c, t) => (CNOT.matrix.clone(), vec![*c, *t]),
+                GateOp::CZ(c, t) => (CZ.matrix.clone(), vec![*c, *t]),
+                GateOp::SWAP(a, b) => (SWAP.matrix.clone(), vec![*a, *b]),
+                GateOp::CCNOT(c1, c2, t) => (TOFFOLI.matrix.clone(), vec![*c1, *c2, *t]),
+                GateOp::CSWAP(c, t1, t2) => (FREDKIN.matrix.clone(), vec![*c, *t1, *t2]),
+                
+                // Non-Clifford fixed gates
+                GateOp::T(t) => (T_GATE.matrix.clone(), vec![*t]),
+                GateOp::Sdg(t) => (SDG_GATE.matrix.clone(), vec![*t]),
+                GateOp::Tdg(t) => (TDG_GATE.matrix.clone(), vec![*t]),
+                GateOp::Sx(t) => (SX_GATE.matrix.clone(), vec![*t]),
+                GateOp::Sxdg(t) => (SXDG_GATE.matrix.clone(), vec![*t]),
+                
+                // Parametric single-qubit gates
+                GateOp::Rx(t, theta) => (rx_matrix(*theta), vec![*t]),
+                GateOp::Ry(t, theta) => (ry_matrix(*theta), vec![*t]),
+                GateOp::Rz(t, theta) => (rz_matrix(*theta), vec![*t]),
+                GateOp::P(t, theta) => (p_matrix(*theta), vec![*t]),
+                GateOp::U1(t, lambda) => (u1_matrix(*lambda), vec![*t]),
+                GateOp::U2(t, phi, lambda) => (u2_matrix(*phi, *lambda), vec![*t]),
+                GateOp::U3(t, theta, phi, lambda) => (u3_matrix(*theta, *phi, *lambda), vec![*t]),
+                
+                // Controlled parametric gates
+                GateOp::CRx(c, t, theta) => (crx_matrix(*theta), vec![*c, *t]),
+                GateOp::CRy(c, t, theta) => (cry_matrix(*theta), vec![*c, *t]),
+                GateOp::CRz(c, t, theta) => (crz_matrix(*theta), vec![*c, *t]),
+                GateOp::CP(c, t, theta) => (cp_matrix(*theta), vec![*c, *t]),
+                
+                // Measurement (skip) and custom gates
                 GateOp::Measure(_, _) => continue,
                 GateOp::Custom(custom_gate, tgts) => {
                     let quantum_gate = custom_gate.to_quantum_gate();
@@ -103,7 +233,7 @@ impl Runtime {
                 }
             };
 
-            state = apply_gate_parallel(&state, &gate.matrix, &targets, num_qubits);
+            state = apply_gate_parallel(&state, &gate_matrix, &targets, num_qubits);
         }
 
         QuantumState::new(state)