[add]: complete gates and functional circuits

5 files changed, 569 insertions, 31 deletions

diff --git a/libpsi-core/src/core/circuit.rs b/libpsi-core/src/core/circuit.rs
index 27c106d..293d509 100644
--- a/libpsi-core/src/core/circuit.rs
+++ b/libpsi-core/src/core/circuit.rs
@@ -1,19 +1,219 @@
-use super::{ClassicalRegister, QuantumRegister};
+use super::{QuantumGate, QuantumRegister, QuantumState};
+use crate::Vector;
+use core::fmt;
+
+#[derive(Clone)]
+pub struct CircuitOperation<'a> {
+    pub gate: &'a QuantumGate<'a>,
+    pub targets: Vec<usize>,
+}
+
+impl<'a> CircuitOperation<'a> {
+    pub fn new(gate: &'a QuantumGate<'a>, targets: Vec<usize>) -> Self {
+        CircuitOperation { gate, targets }
+    }
+}
 
-#[allow(unused)]
 pub struct QuantumCircuit<'a> {
-    quantum_registers: &'a [QuantumRegister<'a>],
-    classical_registers: &'a [ClassicalRegister<'a>],
+    register: QuantumRegister<'a>,
+    operations: Vec<CircuitOperation<'a>>,
 }
 
 impl<'a> QuantumCircuit<'a> {
-    pub fn new(
-        quantum_registers: &'a [QuantumRegister<'a>],
-        classical_registers: &'a [ClassicalRegister<'a>],
-    ) -> QuantumCircuit<'a> {
+    pub fn new(num_qubits: usize) -> QuantumCircuit<'a> {
+        let names: Vec<String> = (0..num_qubits).map(|i| format!("q{}", i)).collect();
+        let leaked_names: &'a [String] = Box::leak(names.into_boxed_slice());
+        let name_refs: Vec<&'a str> = leaked_names.iter().map(|s| s.as_str()).collect();
+
         QuantumCircuit {
-            quantum_registers,
-            classical_registers,
+            register: QuantumRegister::new(
+                Box::leak(Box::new("circuit".to_string())).as_str(),
+                &name_refs,
+            ),
+            operations: Vec::new(),
+        }
+    }
+
+    pub fn from_register(register: QuantumRegister<'a>) -> QuantumCircuit<'a> {
+        QuantumCircuit {
+            register,
+            operations: Vec::new(),
+        }
+    }
+
+    pub fn num_qubits(&self) -> usize {
+        self.register.num_qubits()
+    }
+
+    pub fn state(&self) -> QuantumState {
+        self.register.get_state()
+    }
+
+    pub fn register(&self) -> &QuantumRegister<'a> {
+        &self.register
+    }
+
+    pub fn operations(&self) -> &[CircuitOperation<'a>] {
+        &self.operations
+    }
+
+    pub fn apply(&mut self, gate: &'a QuantumGate<'a>, targets: &[usize]) -> &mut Self {
+        self.register.apply_gate(gate, targets);
+        self.operations
+            .push(CircuitOperation::new(gate, targets.to_vec()));
+        self
+    }
+
+    pub fn h(&mut self, target: usize) -> &mut Self {
+        use crate::gates::HADAMARD;
+        self.register.apply_gate(&HADAMARD, &[target]);
+        self.operations
+            .push(CircuitOperation::new(&HADAMARD, vec![target]));
+        self
+    }
+
+    pub fn x(&mut self, target: usize) -> &mut Self {
+        use crate::gates::PAULI_X;
+        self.register.apply_gate(&PAULI_X, &[target]);
+        self.operations
+            .push(CircuitOperation::new(&PAULI_X, vec![target]));
+        self
+    }
+
+    pub fn y(&mut self, target: usize) -> &mut Self {
+        use crate::gates::PAULI_Y;
+        self.register.apply_gate(&PAULI_Y, &[target]);
+        self.operations
+            .push(CircuitOperation::new(&PAULI_Y, vec![target]));
+        self
+    }
+
+    pub fn z(&mut self, target: usize) -> &mut Self {
+        use crate::gates::PAULI_Z;
+        self.register.apply_gate(&PAULI_Z, &[target]);
+        self.operations
+            .push(CircuitOperation::new(&PAULI_Z, vec![target]));
+        self
+    }
+
+    pub fn s(&mut self, target: usize) -> &mut Self {
+        use crate::gates::S_GATE;
+        self.register.apply_gate(&S_GATE, &[target]);
+        self.operations
+            .push(CircuitOperation::new(&S_GATE, vec![target]));
+        self
+    }
+
+    pub fn t(&mut self, target: usize) -> &mut Self {
+        use crate::gates::T_GATE;
+        self.register.apply_gate(&T_GATE, &[target]);
+        self.operations
+            .push(CircuitOperation::new(&T_GATE, vec![target]));
+        self
+    }
+
+    pub fn cnot(&mut self, control: usize, target: usize) -> &mut Self {
+        use crate::gates::CNOT;
+        self.register.apply_gate(&CNOT, &[control, target]);
+        self.operations
+            .push(CircuitOperation::new(&CNOT, vec![control, target]));
+        self
+    }
+
+    pub fn cx(&mut self, control: usize, target: usize) -> &mut Self {
+        self.cnot(control, target)
+    }
+
+    pub fn cz(&mut self, control: usize, target: usize) -> &mut Self {
+        use crate::gates::CZ;
+        self.register.apply_gate(&CZ, &[control, target]);
+        self.operations
+            .push(CircuitOperation::new(&CZ, vec![control, target]));
+        self
+    }
+
+    pub fn swap(&mut self, qubit1: usize, qubit2: usize) -> &mut Self {
+        use crate::gates::SWAP;
+        self.register.apply_gate(&SWAP, &[qubit1, qubit2]);
+        self.operations
+            .push(CircuitOperation::new(&SWAP, vec![qubit1, qubit2]));
+        self
+    }
+
+    pub fn ccnot(&mut self, control1: usize, control2: usize, target: usize) -> &mut Self {
+        use crate::gates::TOFFOLI;
+        self.register
+            .apply_gate(&TOFFOLI, &[control1, control2, target]);
+        self.operations.push(CircuitOperation::new(
+            &TOFFOLI,
+            vec![control1, control2, target],
+        ));
+        self
+    }
+
+    pub fn toffoli(&mut self, control1: usize, control2: usize, target: usize) -> &mut Self {
+        self.ccnot(control1, control2, target)
+    }
+
+    pub fn cswap(&mut self, control: usize, target1: usize, target2: usize) -> &mut Self {
+        use crate::gates::FREDKIN;
+        self.register
+            .apply_gate(&FREDKIN, &[control, target1, target2]);
+        self.operations.push(CircuitOperation::new(
+            &FREDKIN,
+            vec![control, target1, target2],
+        ));
+        self
+    }
+
+    pub fn fredkin(&mut self, control: usize, target1: usize, target2: usize) -> &mut Self {
+        self.cswap(control, target1, target2)
+    }
+
+    pub fn reset(&mut self) -> &mut Self {
+        let n = self.num_qubits();
+        let names: Vec<String> = (0..n).map(|i| format!("q{}", i)).collect();
+        let leaked_names: &'a [String] = Box::leak(names.into_boxed_slice());
+        let name_refs: Vec<&'a str> = leaked_names.iter().map(|s| s.as_str()).collect();
+
+        self.register = QuantumRegister::new(
+            Box::leak(Box::new("circuit".to_string())).as_str(),
+            &name_refs,
+        );
+        self.operations.clear();
+        self
+    }
+
+    pub fn probability(&self, state_index: usize) -> f64 {
+        let state = self.state();
+        let amp = state.get(state_index);
+        amp.norm2()
+    }
+
+    pub fn probabilities(&self) -> Vec<f64> {
+        let state = self.state();
+        let n = 1 << self.num_qubits();
+        (0..n).map(|i| state.get(i).norm2()).collect()
+    }
+}
+
+impl<'a> fmt::Display for QuantumCircuit<'a> {
+    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
+        writeln!(f, "QuantumCircuit ({} qubits)", self.num_qubits())?;
+        writeln!(f, "Operations:")?;
+        for (i, op) in self.operations.iter().enumerate() {
+            writeln!(f, "  {}: {} on {:?}", i, op.gate.name, op.targets)?;
+        }
+        writeln!(f, "State:")?;
+        let state = self.state();
+        let n = 1 << self.num_qubits();
+        for i in 0..n {
+            let amp = state.get(i);
+            if amp.real.abs() > 1e-10 || amp.imaginary.abs() > 1e-10 {
+                let basis: String = format!("{:0width$b}", i, width = self.num_qubits());
+                writeln!(f, "  |{}⟩: {:.4}", basis, amp)?;
+            }
         }
+        Ok(())
     }
 }
diff --git a/libpsi-core/src/core/classical_components.rs b/libpsi-core/src/core/classical_components.rs
index 6710e9d..f6565af 100644
--- a/libpsi-core/src/core/classical_components.rs
+++ b/libpsi-core/src/core/classical_components.rs
@@ -13,7 +13,7 @@ pub struct ClassicalRegister<'a> {
 }
 
 impl<'a> ClassicalBit<'a> {
-    pub fn new(name: &'a str, state: bool) -> ClassicalBit {
+    pub fn new(name: &'a str, state: bool) -> ClassicalBit<'a> {
         ClassicalBit { name, state }
     }
 
diff --git a/libpsi-core/src/core/gates.rs b/libpsi-core/src/core/gates.rs
index 854313d..4f922a8 100644
--- a/libpsi-core/src/core/gates.rs
+++ b/libpsi-core/src/core/gates.rs
@@ -6,25 +6,50 @@ lazy_static::lazy_static! {
         name: "H",
         matrix: matrix!([complex!(1.0, 0.0), complex!( 1.0, 0.0)];
                         [complex!(1.0, 0.0), complex!(-1.0, 0.0)]) *
-                complex!(1.0/2.0_f64.sqrt(), 0.0)
+                complex!(1.0/2.0_f64.sqrt(), 0.0),
+        num_qubits: 1,
     };
 
     pub static ref PAULI_X: QuantumGate<'static> = QuantumGate {
-        name: "Pauli-X",
+        name: "X",
         matrix: matrix!([complex!(0.0, 0.0), complex!(1.0, 0.0)];
-                        [complex!(1.0, 0.0), complex!(0.0, 0.0)])
+                        [complex!(1.0, 0.0), complex!(0.0, 0.0)]),
+        num_qubits: 1,
     };
 
     pub static ref PAULI_Y: QuantumGate<'static> = QuantumGate {
-        name: "Pauli-Y", 
+        name: "Y", 
         matrix: matrix!([complex!(0.0, 0.0), complex!(0.0, -1.0)];
-                        [complex!(0.0, 1.0), complex!(0.0,  0.0)])
+                        [complex!(0.0, 1.0), complex!(0.0,  0.0)]),
+        num_qubits: 1,
     };
 
     pub static ref PAULI_Z: QuantumGate<'static> = QuantumGate {
-        name: "Pauli-Z", 
+        name: "Z", 
         matrix: matrix!([complex!(1.0, 0.0), complex!( 0.0, 0.0)];
-                        [complex!(0.0, 0.0), complex!(-1.0, 0.0)])
+                        [complex!(0.0, 0.0), complex!(-1.0, 0.0)]),
+        num_qubits: 1,
+    };
+    
+    pub static ref S_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 T_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 IDENTITY: QuantumGate<'static> = QuantumGate {
+        name: "I",
+        matrix: matrix!([complex!(1.0, 0.0), complex!(0.0, 0.0)];
+                        [complex!(0.0, 0.0), complex!(1.0, 0.0)]),
+        num_qubits: 1,
     };
 
     pub static ref CNOT: QuantumGate<'static> = QuantumGate {
@@ -32,6 +57,73 @@ lazy_static::lazy_static! {
         matrix: 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!(0.0, 0.0), complex!(1.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!(1.0, 0.0), complex!(0.0, 0.0)]),
+        num_qubits: 2,
+    };
+    
+    pub static ref CZ: QuantumGate<'static> = QuantumGate {
+        name: "CZ",
+        matrix: 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!(-1.0, 0.0)]),
+        num_qubits: 2,
+    };
+    
+    pub static ref SWAP: QuantumGate<'static> = QuantumGate {
+        name: "SWAP",
+        matrix: 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!(0.0, 0.0), complex!(1.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!(0.0, 0.0), complex!(1.0, 0.0)]),
+        num_qubits: 2,
+    };
+    
+    pub static ref ISWAP: QuantumGate<'static> = QuantumGate {
+        name: "iSWAP",
+        matrix: 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!(0.0, 0.0), complex!(0.0, 1.0), complex!(0.0, 0.0)];
+                        [complex!(0.0, 0.0), complex!(0.0, 1.0), complex!(0.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)]),
+        num_qubits: 2,
+    };
+    
+    pub static ref SQRT_SWAP: QuantumGate<'static> = QuantumGate {
+        name: "√SWAP",
+        matrix: 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!(0.5, 0.5), complex!(0.5, -0.5), complex!(0.0, 0.0)];
+                        [complex!(0.0, 0.0), complex!(0.5, -0.5), complex!(0.5, 0.5), 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)]),
+        num_qubits: 2,
+    };
+
+    pub static ref TOFFOLI: QuantumGate<'static> = QuantumGate {
+        name: "CCNOT",
+        matrix: 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!(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!(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!(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!(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!(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!(0.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!(0.0, 0.0), complex!(0.0, 0.0), complex!(1.0, 0.0), complex!(0.0, 0.0)]
+        ),
+        num_qubits: 3,
+    };
+    
+    pub static ref FREDKIN: QuantumGate<'static> = QuantumGate {
+        name: "CSWAP",
+        matrix: 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!(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!(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!(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!(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!(0.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!(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!(0.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)]
+        ),
+        num_qubits: 3,
     };
 }
diff --git a/libpsi-core/src/core/quantum_components.rs b/libpsi-core/src/core/quantum_components.rs
index e6931a0..75f6325 100644
--- a/libpsi-core/src/core/quantum_components.rs
+++ b/libpsi-core/src/core/quantum_components.rs
@@ -1,7 +1,6 @@
-use crate::{column_vector, complex, ColumnVector, Complex, Matrix, Vector, VectorMatrix};
+use crate::{column_vector, complex, ColumnVector, Complex, Float, Matrix, Vector, VectorMatrix};
 use core::{fmt, ops};
 
-// TODO(Hachem): Redo these macros to work with the new function definition.
 #[macro_export]
 macro_rules! count {
     () => { 0 };
@@ -43,6 +42,14 @@ impl QuantumState {
     }
 }
 
+fn identity_matrix<T: Float>(size: usize) -> Matrix<T> {
+    let mut data = vec![T::zero(); size * size];
+    for i in 0..size {
+        data[i * size + i] = T::one();
+    }
+    Matrix::new(size, size, data)
+}
+
 #[derive(Clone)]
 pub struct QuantumBit<'a> {
     state: QuantumState,
@@ -60,6 +67,41 @@ pub struct QuantumRegister<'a> {
 pub struct QuantumGate<'a> {
     pub name: &'a str,
     pub matrix: Matrix<Complex<f64>>,
+    pub num_qubits: usize,
+}
+
+impl<'a> QuantumGate<'a> {
+    pub fn new(name: &'a str, matrix: Matrix<Complex<f64>>, num_qubits: usize) -> Self {
+        let expected_dim = 1 << num_qubits;
+        assert_eq!(
+            matrix.rows, expected_dim,
+            "Gate matrix rows must be 2^num_qubits"
+        );
+        assert_eq!(
+            matrix.cols, expected_dim,
+            "Gate matrix cols must be 2^num_qubits"
+        );
+        QuantumGate {
+            name,
+            matrix,
+            num_qubits,
+        }
+    }
+
+    pub fn from_matrix(name: &'a str, matrix: Matrix<Complex<f64>>) -> Self {
+        assert_eq!(matrix.rows, matrix.cols, "Gate matrix must be square");
+        let dim = matrix.rows;
+        assert!(
+            dim > 0 && (dim & (dim - 1)) == 0,
+            "Matrix dimension must be a power of 2"
+        );
+        let num_qubits = (dim as f64).log2() as usize;
+        QuantumGate {
+            name,
+            matrix,
+            num_qubits,
+        }
+    }
 }
 
 impl<'a> QuantumBit<'a> {
@@ -111,7 +153,7 @@ impl<'a> QuantumRegister<'a> {
         self.state_vector = ColumnVector::from_matrix(&new_result);
     }
 
-    pub fn get_bits(&self) -> Vec<QuantumBit> {
+    pub fn get_bits(&self) -> Vec<QuantumBit<'_>> {
         self.qubits.clone()
     }
 
@@ -122,6 +164,137 @@ impl<'a> QuantumRegister<'a> {
     pub fn get_name(&self) -> &'a str {
         self.name
     }
+
+    pub fn num_qubits(&self) -> usize {
+        self.qubits.len()
+    }
+
+    pub fn apply_gate(&mut self, gate: &QuantumGate, targets: &[usize]) {
+        let n = self.num_qubits();
+
+        assert_eq!(
+            gate.num_qubits,
+            targets.len(),
+            "Number of target qubits must match gate's qubit count"
+        );
+        for &t in targets {
+            assert!(
+                t < n,
+                "Target qubit index {} out of range for {}-qubit register",
+                t,
+                n
+            );
+        }
+
+        let mut sorted_targets = targets.to_vec();
+        sorted_targets.sort();
+        for i in 1..sorted_targets.len() {
+            assert_ne!(
+                sorted_targets[i],
+                sorted_targets[i - 1],
+                "Duplicate target qubit indices are not allowed"
+            );
+        }
+
+        let full_operator = self.build_full_operator(gate, targets);
+
+        self.state_vector = self
+            .state_vector
+            .mul_matrix(&full_operator)
+            .expect("Matrix multiplication failed during gate application");
+    }
+
+    fn build_full_operator(&self, gate: &QuantumGate, targets: &[usize]) -> Matrix<Complex<f64>> {
+        let n = self.num_qubits();
+        let g = gate.num_qubits;
+        let dim = 1 << n;
+
+        let mut contiguous = true;
+        for i in 1..targets.len() {
+            if targets[i] != targets[i - 1] + 1 {
+                contiguous = false;
+                break;
+            }
+        }
+
+        if contiguous && g == n {
+            return gate.matrix.clone();
+        }
+
+        if contiguous {
+            return self.build_contiguous_operator(gate, targets[0]);
+        }
+
+        let mut result = Matrix::new(dim, dim, vec![complex!(0.0, 0.0); dim * dim]);
+
+        for col in 0..dim {
+            for row in 0..dim {
+                let mut target_row_bits = 0usize;
+                let mut target_col_bits = 0usize;
+
+                for (i, &t) in targets.iter().enumerate() {
+                    let qubit_pos = n - 1 - t;
+                    if (row >> qubit_pos) & 1 == 1 {
+                        target_row_bits |= 1 << (g - 1 - i);
+                    }
+                    if (col >> qubit_pos) & 1 == 1 {
+                        target_col_bits |= 1 << (g - 1 - i);
+                    }
+                }
+
+                let mut non_target_match = true;
+                for q in 0..n {
+                    if !targets.contains(&q) {
+                        let qubit_pos = n - 1 - q;
+                        if ((row >> qubit_pos) & 1) != ((col >> qubit_pos) & 1) {
+                            non_target_match = false;
+                            break;
+                        }
+                    }
+                }
+
+                if non_target_match {
+                    result.set(row, col, gate.matrix.get(target_row_bits, target_col_bits));
+                }
+            }
+        }
+
+        result
+    }
+
+    fn build_contiguous_operator(
+        &self,
+        gate: &QuantumGate,
+        start_idx: usize,
+    ) -> Matrix<Complex<f64>> {
+        let n = self.num_qubits();
+        let g = gate.num_qubits;
+
+        let mut result: Option<Matrix<Complex<f64>>> = None;
+
+        for i in 0..n {
+            let part: Matrix<Complex<f64>> = if i == start_idx {
+                gate.matrix.clone()
+            } else if i > start_idx && i < start_idx + g {
+                continue;
+            } else {
+                identity_matrix(2)
+            };
+
+            result = Some(match result {
+                None => part,
+                Some(r) => r.kronecker(&part),
+            });
+        }
+
+        result.unwrap_or_else(|| identity_matrix(1 << n))
+    }
+
+    pub fn apply_gates(&mut self, operations: &[(&QuantumGate, &[usize])]) {
+        for (gate, targets) in operations {
+            self.apply_gate(gate, targets);
+        }
+    }
 }
 
 impl<'a> ops::Index<usize> for QuantumRegister<'a> {
diff --git a/tester/src/main.rs b/tester/src/main.rs
index 2c2e72c..96fc86f 100644
--- a/tester/src/main.rs
+++ b/tester/src/main.rs
@@ -1,15 +1,88 @@
 use libpsi_core::*;
 
 fn main() {
-    let quantum_registers = [
-        QuantumRegister::new("qr0", &["q0", "q1", "q2", "q3"]),
-        QuantumRegister::new("qr1", &["k0", "k1", "k2", "k3"]),
-    ];
+    println!("Bell State Creation: |Φ+⟩ = (|00⟩ + |11⟩)/√2");
+    println!("   Circuit: H(q0) → CNOT(q0, q1)");
+    let mut bell = QuantumCircuit::new(2);
+    bell.h(0).cnot(0, 1);
+    println!("{}", bell);
 
-    let classical_registers = [
-        ClassicalRegister::new("cr0", &["c0", "c1", "c2", "c3"]),
-        ClassicalRegister::new("cr1", &["a0", "a1", "a2", "a3"]),
-    ];
+    print!("\n------\n\n");
 
-    QuantumCircuit::new(&quantum_registers, &classical_registers);
+    println!("GHZ State (3-qubit entanglement): |GHZ> = (|000⟩ + |111⟩)/√2");
+    println!("   Circuit: H(q0) → CNOT(q0, q1) → CNOT(q0, q2)");
+    let mut ghz = QuantumCircuit::new(3);
+    ghz.h(0).cnot(0, 1).cnot(0, 2);
+    println!("{}", ghz);
+
+    print!("\n------\n\n");
+
+    println!("SWAP via 3 CNOTs");
+    println!("   Start with |10>, apply CNOT chain");
+    let mut swap_circuit = QuantumCircuit::new(2);
+    swap_circuit
+        .x(0) // Set to |10⟩
+        .cnot(0, 1)
+        .cnot(1, 0)
+        .cnot(0, 1);
+    println!("{}", swap_circuit);
+
+    print!("\n------\n\n");
+
+    println!("Toffoli Gate (Reversible AND)");
+    println!("   CCNOT flips q2 only when q0=1 AND q1=1");
+    let mut toffoli_circuit = QuantumCircuit::new(3);
+    toffoli_circuit
+        .x(0)
+        .x(1) // Set to |110⟩
+        .toffoli(0, 1, 2);
+    println!("{}", toffoli_circuit);
+
+    print!("\n------\n\n");
+
+    println!("Fredkin Gate (Controlled SWAP)");
+    println!("   CSWAP swaps q1 and q2 only when q0=1");
+    let mut fredkin_circuit = QuantumCircuit::new(3);
+    fredkin_circuit
+        .x(0)
+        .x(1) // Set to |110⟩
+        .fredkin(0, 1, 2);
+    println!("{}", fredkin_circuit);
+
+    print!("\n------\n\n");
+
+    println!("6. Non-contiguous CNOT (q0 controls q2, skipping q1)");
+    let mut nc_circuit = QuantumCircuit::new(3);
+    nc_circuit
+        .x(0) // |100⟩
+        .cnot(0, 2); // CNOT with control=q0, target=q2
+    println!("{}", nc_circuit);
+
+    print!("\n------\n\n");
+
+    println!("Full Superposition (H on all qubits)");
+    let mut super_circuit = QuantumCircuit::new(3);
+    super_circuit.h(0).h(1).h(2);
+    println!("{}", super_circuit);
+
+    print!("\n------\n\n");
+
+    println!("Probability Test");
+    let mut prob_circuit = QuantumCircuit::new(2);
+    prob_circuit.h(0).cnot(0, 1);
+    println!("   Bell state probabilities:");
+    let probs = prob_circuit.probabilities();
+    for (i, p) in probs.iter().enumerate() {
+        if *p > 1e-10 {
+            println!("     |{:02b}⟩: {:.4}", i, p);
+        }
+    }
+    println!();
+
+    print!("\n------\n\n");
+
+    println!("Complex Circuit with Method Chaining");
+    let mut complex = QuantumCircuit::new(4);
+    complex.h(0).h(1).cnot(0, 2).cnot(1, 3).cz(2, 3).swap(0, 1);
+    println!("{}", complex);
 }