feat: simd + testing + formatting

18 files changed, 449 insertions, 1829 deletions

diff --git a/tester/benchmarks.rs b/tester/benchmarks.rs
deleted file mode 100644
index 7aea535..0000000
--- a/tester/benchmarks.rs
+++ /dev/null
@@ -1,90 +0,0 @@
-use crate::common::{benchmark_circuit, print_section, BenchmarkResult};
-use psi::QuantumCircuit;
-use psi::HorizontalRenderer;
-
-pub fn run_all(results: &mut Vec<BenchmarkResult>) {
-    println!("═══════════════════════════════════════════════════════════════");
-    println!("                    BENCHMARK CIRCUITS");
-    println!("═══════════════════════════════════════════════════════════════\n");
-
-    test_8_qubit(results);
-    test_10_qubit(results);
-    test_12_qubit(results);
-    test_14_qubit(results);
-}
-
-pub fn test_8_qubit(results: &mut Vec<BenchmarkResult>) {
-    print_section("8-qubit Entangled Circuit");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(8);
-        for i in 0..8 {
-            circuit.h(i);
-        }
-        for i in 0..7 {
-            circuit.cnot(i, i + 1);
-        }
-        circuit
-    };
-
-    println!("{}", HorizontalRenderer::new(&builder()));
-    results.push(benchmark_circuit("8-qubit entangled", builder));
-}
-
-pub fn test_10_qubit(results: &mut Vec<BenchmarkResult>) {
-    print_section("10-qubit Entangled Circuit");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(10);
-        for i in 0..10 {
-            circuit.h(i);
-        }
-        for i in 0..9 {
-            circuit.cnot(i, i + 1);
-        }
-        circuit.cz(0, 9);
-        circuit
-    };
-
-    println!("{}", HorizontalRenderer::new(&builder()));
-    results.push(benchmark_circuit("10-qubit entangled", builder));
-}
-
-pub fn test_12_qubit(results: &mut Vec<BenchmarkResult>) {
-    print_section("12-qubit Entangled Circuit");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(12);
-        for i in 0..12 {
-            circuit.h(i);
-        }
-        for i in 0..11 {
-            circuit.cnot(i, i + 1);
-        }
-        circuit.cz(0, 11);
-        circuit.swap(5, 6);
-        circuit
-    };
-
-    println!("{}", HorizontalRenderer::new(&builder()));
-    results.push(benchmark_circuit("12-qubit entangled", builder));
-}
-
-pub fn test_14_qubit(results: &mut Vec<BenchmarkResult>) {
-    print_section("14-qubit Entangled Circuit");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(14);
-        for i in 0..14 {
-            circuit.h(i);
-        }
-        for i in 0..13 {
-            circuit.cnot(i, i + 1);
-        }
-        circuit
-    };
-
-    println!("{}", HorizontalRenderer::new(&builder()));
-    results.push(benchmark_circuit("14-qubit entangled", builder));
-}
-
diff --git a/tester/clifford.c b/tester/clifford.c
new file mode 100644
index 0000000..3be9234
--- /dev/null
+++ b/tester/clifford.c
@@ -0,0 +1,75 @@
+#include "tests.h"
+
+static const double R2 = 0.7071067811865476;
+
+static void build_bell(struct PsiQuantumCircuit* c)
+{
+	psi_apply_h(c, 0);
+	psi_apply_cnot(c, 0, 1);
+}
+
+static void build_ghz(struct PsiQuantumCircuit* c)
+{
+	psi_apply_h(c, 0);
+	psi_apply_cnot(c, 0, 1);
+	psi_apply_cnot(c, 0, 2);
+}
+
+void run_clifford_tests(void)
+{
+	psi_test_section("Clifford gates");
+
+	struct PsiQuantumCircuit x = psi_new_quantum_circuit(1);
+	psi_apply_x(&x, 0);
+	struct PsiComplex x_exp[] = { psi_new_complex(0.0, 0.0), psi_new_complex(1.0, 0.0) };
+	psi_check_circuit("X|0> = |1>", &x, x_exp, 2);
+	psi_free_quantum_circuit(&x);
+
+	struct PsiQuantumCircuit bell = psi_new_quantum_circuit(2);
+	build_bell(&bell);
+	struct PsiComplex bell_exp[] = {
+		psi_new_complex(R2, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(R2, 0.0),
+	};
+	psi_check_circuit("Bell = (|00>+|11>)/sqrt2", &bell, bell_exp, 4);
+	psi_free_quantum_circuit(&bell);
+
+	struct PsiQuantumCircuit ghz = psi_new_quantum_circuit(3);
+	build_ghz(&ghz);
+	struct PsiComplex ghz_exp[] = {
+		psi_new_complex(R2, 0.0),  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		psi_new_complex(0.0, 0.0), psi_new_complex(R2, 0.0),
+	};
+	psi_check_circuit("GHZ = (|000>+|111>)/sqrt2", &ghz, ghz_exp, 8);
+	psi_free_quantum_circuit(&ghz);
+
+	struct PsiQuantumCircuit swap = psi_new_quantum_circuit(2);
+	psi_apply_x(&swap, 0);
+	psi_apply_cnot(&swap, 0, 1);
+	psi_apply_cnot(&swap, 1, 0);
+	psi_apply_cnot(&swap, 0, 1);
+	struct PsiComplex swap_exp[] = {
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(1.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+	};
+	psi_check_circuit("SWAP via 3 CNOTs: |10> -> |01>", &swap, swap_exp, 4);
+	psi_free_quantum_circuit(&swap);
+
+	struct PsiQuantumCircuit toffoli = psi_new_quantum_circuit(3);
+	psi_apply_x(&toffoli, 0);
+	psi_apply_x(&toffoli, 1);
+	psi_apply_ccnot(&toffoli, 0, 1, 2);
+	struct PsiComplex toffoli_exp[8];
+	for (size_t i = 0; i < 8; i++)
+		toffoli_exp[i] = psi_new_complex(0.0, 0.0);
+	toffoli_exp[7] = psi_new_complex(1.0, 0.0);
+	psi_check_circuit("Toffoli: |110> -> |111>", &toffoli, toffoli_exp, 8);
+	psi_free_quantum_circuit(&toffoli);
+
+	psi_check_runtimes_agree("Runtimes agree on GHZ", 3, build_ghz);
+}
diff --git a/tester/clifford.rs b/tester/clifford.rs
deleted file mode 100644
index 932a918..0000000
--- a/tester/clifford.rs
+++ /dev/null
@@ -1,126 +0,0 @@
-use crate::common::{benchmark_circuit, print_circuit, print_section, BenchmarkResult};
-use psi::QuantumCircuit;
-
-pub fn run_all(results: &mut Vec<BenchmarkResult>) {
-    println!("═══════════════════════════════════════════════════════════════");
-    println!("                     CLIFFORD GATES TESTS");
-    println!("═══════════════════════════════════════════════════════════════\n");
-
-    test_bell_state(results);
-    test_ghz_state(results);
-    test_swap_via_cnots(results);
-    test_toffoli(results);
-    test_hadamard_measure(results);
-    test_complex_circuit(results);
-}
-
-pub fn test_bell_state(results: &mut Vec<BenchmarkResult>) {
-    print_section("Bell State with Measurement");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::with_classical(2, 2);
-        circuit.h(0).cnot(0, 1).measure(0, 0).measure(1, 1);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("Bell State (2 qubits)", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_ghz_state(results: &mut Vec<BenchmarkResult>) {
-    print_section("GHZ State");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(3);
-        circuit.h(0).cnot(0, 1).cnot(0, 2);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("GHZ State (3 qubits)", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_swap_via_cnots(results: &mut Vec<BenchmarkResult>) {
-    print_section("SWAP via 3 CNOTs");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(2);
-        circuit.x(0).cnot(0, 1).cnot(1, 0).cnot(0, 1);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("SWAP via CNOTs (2 qubits)", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_toffoli(results: &mut Vec<BenchmarkResult>) {
-    print_section("Toffoli Gate");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(3);
-        circuit.x(0).x(1).toffoli(0, 1, 2);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("Toffoli (3 qubits)", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_hadamard_measure(results: &mut Vec<BenchmarkResult>) {
-    print_section("Full Circuit with Measurements");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::with_classical(3, 3);
-        circuit.h(0).h(1).h(2).measure_all();
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("3-qubit Hadamard + Measure", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_complex_circuit(results: &mut Vec<BenchmarkResult>) {
-    print_section("Complex Circuit");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::with_classical(4, 2);
-        circuit
-            .h(0)
-            .h(1)
-            .cnot(0, 2)
-            .cnot(1, 3)
-            .cz(2, 3)
-            .swap(0, 1)
-            .measure(0, 0)
-            .measure(1, 1);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("Complex (4 qubits)", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
diff --git a/tester/common.rs b/tester/common.rs
deleted file mode 100644
index 9e576da..0000000
--- a/tester/common.rs
+++ /dev/null
@@ -1,215 +0,0 @@
-use psi::{QuantumCircuit, QuantumState, Runtime, Vector};
-use psi::{HorizontalRenderer, VerticalRenderer};
-use std::time::{Duration, Instant};
-
-/// A named list of circuit builders used by the benchmark/test suites.
-pub type CircuitCases = Vec<(&'static str, Box<dyn Fn() -> QuantumCircuit>)>;
-
-pub struct BenchmarkResult {
-    pub name: String,
-    pub basic_time: Duration,
-    pub mt_time: Duration,
-    pub results_match: bool,
-}
-
-pub fn benchmark_circuit<F>(name: &str, circuit_builder: F) -> BenchmarkResult
-where
-    F: Fn() -> QuantumCircuit,
-{
-    let mut circuit_st = circuit_builder();
-    let mut circuit_mt = circuit_builder();
-
-    let start_st = Instant::now();
-    circuit_st.compute_with(Runtime::BasicRT);
-    let basic_time = start_st.elapsed();
-
-    let start_mt = Instant::now();
-    circuit_mt.compute_with(Runtime::BasicRTMT);
-    let mt_time = start_mt.elapsed();
-
-    let state_st = circuit_st.state();
-    let state_mt = circuit_mt.state();
-
-    let results_match = states_equal(state_st, state_mt);
-
-    BenchmarkResult {
-        name: name.to_string(),
-        basic_time,
-        mt_time,
-        results_match,
-    }
-}
-
-pub fn states_equal(a: &QuantumState, b: &QuantumState) -> bool {
-    if a.size() != b.size() {
-        return false;
-    }
-    for i in 0..a.size() {
-        let amp_a = a.get(i);
-        let amp_b = b.get(i);
-        let diff_real = (amp_a.real - amp_b.real).abs();
-        let diff_imag = (amp_a.imaginary - amp_b.imaginary).abs();
-        if diff_real > 1e-10 || diff_imag > 1e-10 {
-            return false;
-        }
-    }
-    true
-}
-
-pub fn format_duration(d: Duration) -> String {
-    if d.as_secs() > 0 {
-        format!("{:.3}s", d.as_secs_f64())
-    } else if d.as_millis() > 0 {
-        format!("{:.3}ms", d.as_secs_f64() * 1000.0)
-    } else {
-        format!("{:.3}μs", d.as_secs_f64() * 1_000_000.0)
-    }
-}
-
-pub fn print_section(title: &str) {
-    let width = 61;
-    let padding = width - title.len() - 2;
-    println!("┌{}┐", "─".repeat(width));
-    println!("│ {}{} │", title, " ".repeat(padding));
-    println!("└{}┘\n", "─".repeat(width));
-}
-
-pub fn print_circuit(circuit: &QuantumCircuit) {
-    println!("Horizontal:\n{}", HorizontalRenderer::new(circuit));
-    println!("Vertical:\n{}", VerticalRenderer::new(circuit));
-}
-
-pub fn print_benchmark_table(results: &[BenchmarkResult]) {
-    if results.is_empty() {
-        return;
-    }
-
-    let headers = ["Circuit", "BasicRT", "BasicRTMT", "Speedup", "Match"];
-
-    let formatted: Vec<(String, String, String, String, String)> = results
-        .iter()
-        .map(|r| {
-            let speedup = r.basic_time.as_secs_f64() / r.mt_time.as_secs_f64();
-            (
-                r.name.clone(),
-                format_duration(r.basic_time),
-                format_duration(r.mt_time),
-                if speedup.is_finite() {
-                    format!("{:.2}x", speedup)
-                } else {
-                    "N/A".to_string()
-                },
-                if r.results_match { "✓" } else { "✗" }.to_string(),
-            )
-        })
-        .collect();
-
-    let c1 = formatted
-        .iter()
-        .map(|r| r.0.len())
-        .max()
-        .unwrap()
-        .max(headers[0].len());
-    let c2 = formatted
-        .iter()
-        .map(|r| r.1.len())
-        .max()
-        .unwrap()
-        .max(headers[1].len());
-    let c3 = formatted
-        .iter()
-        .map(|r| r.2.len())
-        .max()
-        .unwrap()
-        .max(headers[2].len());
-    let c4 = formatted
-        .iter()
-        .map(|r| r.3.len())
-        .max()
-        .unwrap()
-        .max(headers[3].len());
-    let c5 = formatted
-        .iter()
-        .map(|r| r.4.chars().count())
-        .max()
-        .unwrap()
-        .max(headers[4].len());
-
-    let top = format!(
-        "╔{}═{}═{}═{}═{}╗",
-        "═".repeat(c1 + 2),
-        "═".repeat(c2 + 2),
-        "═".repeat(c3 + 2),
-        "═".repeat(c4 + 2),
-        "═".repeat(c5 + 2)
-    );
-    let title_sep = format!(
-        "╠{}╤{}╤{}╤{}╤{}╣",
-        "═".repeat(c1 + 2),
-        "═".repeat(c2 + 2),
-        "═".repeat(c3 + 2),
-        "═".repeat(c4 + 2),
-        "═".repeat(c5 + 2)
-    );
-    let header_sep = format!(
-        "╠{}╪{}╪{}╪{}╪{}╣",
-        "═".repeat(c1 + 2),
-        "═".repeat(c2 + 2),
-        "═".repeat(c3 + 2),
-        "═".repeat(c4 + 2),
-        "═".repeat(c5 + 2)
-    );
-    let bottom = format!(
-        "╚{}╧{}╧{}╧{}╧{}╝",
-        "═".repeat(c1 + 2),
-        "═".repeat(c2 + 2),
-        "═".repeat(c3 + 2),
-        "═".repeat(c4 + 2),
-        "═".repeat(c5 + 2)
-    );
-
-    let total_width = c1 + c2 + c3 + c4 + c5 + 14;
-
-    println!("\n{}", top);
-    println!(
-        "║{:^width$}║",
-        "RUNTIME BENCHMARK RESULTS",
-        width = total_width
-    );
-    println!("{}", title_sep);
-    println!(
-        "║ {:<c1$} │ {:^c2$} │ {:^c3$} │ {:^c4$} │ {:^c5$} ║",
-        headers[0], headers[1], headers[2], headers[3], headers[4],
-    );
-    println!("{}", header_sep);
-
-    for (name, basic, mt, speedup, matched) in &formatted {
-        println!(
-            "║ {:<c1$} │ {:>c2$} │ {:>c3$} │ {:>c4$} │ {:^c5$} ║",
-            name, basic, mt, speedup, matched,
-        );
-    }
-
-    println!("{}", bottom);
-}
-
-pub fn print_summary(results: &[BenchmarkResult]) {
-    let all_match = results.iter().all(|r| r.results_match);
-    println!("\n");
-    if all_match {
-        println!("✓ All circuits produced identical results with both runtimes!");
-    } else {
-        println!("✗ WARNING: Some circuits produced different results!");
-    }
-
-    let total_basic: Duration = results.iter().map(|r| r.basic_time).sum();
-    let total_mt: Duration = results.iter().map(|r| r.mt_time).sum();
-    let overall_speedup = total_basic.as_secs_f64() / total_mt.as_secs_f64();
-
-    println!(
-        "\nTotal time - BasicRT: {} | BasicRTMT: {} | Overall speedup: {:.2}x",
-        format_duration(total_basic),
-        format_duration(total_mt),
-        overall_speedup
-    );
-}
diff --git a/tester/custom_gates.c b/tester/custom_gates.c
new file mode 100644
index 0000000..df2d9ae
--- /dev/null
+++ b/tester/custom_gates.c
@@ -0,0 +1,36 @@
+#include "tests.h"
+
+static const double R2 = 0.7071067811865476;
+
+void run_custom_tests(void)
+{
+	psi_test_section("Custom gates");
+
+	struct PsiCompositeGateOp bell_ops[] = {
+		{ PSI_OP_H, { 0 }, 1 },
+		{ PSI_OP_CNOT, { 0, 1 }, 2 },
+	};
+	struct PsiCustomGate bell = psi_new_custom_gate_from_composite("BELL", 2, bell_ops, 2);
+	struct PsiQuantumCircuit c = psi_new_quantum_circuit(2);
+	size_t bell_targets[] = { 0, 1 };
+	psi_apply_custom(&c, bell, bell_targets, 2);
+	struct PsiComplex bell_exp[] = {
+		psi_new_complex(R2, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(R2, 0.0),
+	};
+	psi_check_circuit("composite BELL gate -> Bell state", &c, bell_exp, 4);
+	psi_free_quantum_circuit(&c);
+
+	struct PsiCustomGate custom_x = psi_new_custom_gate_from_matrix(
+			"MYX",
+			psi_matrix(2, 2, psi_new_complex(0.0, 0.0), psi_new_complex(1.0, 0.0),
+			           psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0)));
+	struct PsiQuantumCircuit c2 = psi_new_quantum_circuit(1);
+	size_t x_targets[] = { 0 };
+	psi_apply_custom(&c2, custom_x, x_targets, 1);
+	struct PsiComplex x_exp[] = { psi_new_complex(0.0, 0.0), psi_new_complex(1.0, 0.0) };
+	psi_check_circuit("from_matrix X custom gate: |0> -> |1>", &c2, x_exp, 2);
+	psi_free_quantum_circuit(&c2);
+}
diff --git a/tester/custom_gates.rs b/tester/custom_gates.rs
deleted file mode 100644
index ffe5c01..0000000
--- a/tester/custom_gates.rs
+++ /dev/null
@@ -1,121 +0,0 @@
-use crate::common::{benchmark_circuit, print_circuit, print_section, BenchmarkResult};
-use psi::{complex, matrix, CustomGate, CustomGateBuilder, QuantumCircuit};
-
-pub fn run_all(results: &mut Vec<BenchmarkResult>) {
-    println!("═══════════════════════════════════════════════════════════════");
-    println!("                    CUSTOM GATES TESTS");
-    println!("═══════════════════════════════════════════════════════════════\n");
-
-    test_bell_gate(results);
-    test_swap_gate(results);
-    test_sqrt_x_gate(results);
-}
-
-pub fn test_bell_gate(results: &mut Vec<BenchmarkResult>) {
-    print_section("Custom Gate: Bell Pair Creator");
-
-    let bell_gate = CustomGateBuilder::new("BELL", 2).h(0).cnot(0, 1).build();
-    let gate_clone = bell_gate.clone();
-
-    let builder = move || {
-        let mut circuit = QuantumCircuit::new(4);
-        circuit
-            .apply_custom(gate_clone.clone(), &[0, 1])
-            .apply_custom(gate_clone.clone(), &[2, 3]);
-        circuit
-    };
-
-    let display_circuit = {
-        let mut circuit = QuantumCircuit::new(4);
-        circuit
-            .apply_custom(bell_gate.clone(), &[0, 1])
-            .apply_custom(bell_gate.clone(), &[2, 3]);
-        circuit
-    };
-    print_circuit(&display_circuit);
-    results.push(benchmark_circuit("Custom BELL (4 qubits)", builder));
-
-    let mut display = {
-        let mut circuit = QuantumCircuit::new(4);
-        circuit
-            .apply_custom(bell_gate.clone(), &[0, 1])
-            .apply_custom(bell_gate.clone(), &[2, 3]);
-        circuit
-    };
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_swap_gate(results: &mut Vec<BenchmarkResult>) {
-    print_section("Custom Gate: Swap via CNOTs");
-
-    let swap_gate = CustomGateBuilder::new("MYSWAP", 2)
-        .cnot(0, 1)
-        .cnot(1, 0)
-        .cnot(0, 1)
-        .build();
-    let gate_clone = swap_gate.clone();
-
-    let builder = move || {
-        let mut circuit = QuantumCircuit::new(2);
-        circuit.x(0).apply_custom(gate_clone.clone(), &[0, 1]);
-        circuit
-    };
-
-    let display_circuit = {
-        let mut circuit = QuantumCircuit::new(2);
-        circuit.x(0).apply_custom(swap_gate.clone(), &[0, 1]);
-        circuit
-    };
-    print_circuit(&display_circuit);
-    results.push(benchmark_circuit("Custom SWAP (2 qubits)", builder));
-
-    let mut display = {
-        let mut circuit = QuantumCircuit::new(2);
-        circuit.x(0).apply_custom(swap_gate.clone(), &[0, 1]);
-        circuit
-    };
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_sqrt_x_gate(results: &mut Vec<BenchmarkResult>) {
-    print_section("Custom Gate: Matrix-defined √X gate");
-
-    let sqrt_x_matrix = matrix!(
-        [complex!(0.5, 0.5), complex!(0.5, -0.5)];
-        [complex!(0.5, -0.5), complex!(0.5, 0.5)]
-    );
-    let sqrt_x = CustomGate::from_matrix("√X", sqrt_x_matrix);
-    let gate_clone = sqrt_x.clone();
-
-    let builder = move || {
-        let mut circuit = QuantumCircuit::new(1);
-        circuit
-            .apply_custom(gate_clone.clone(), &[0])
-            .apply_custom(gate_clone.clone(), &[0]);
-        circuit
-    };
-
-    let display_circuit = {
-        let mut circuit = QuantumCircuit::new(1);
-        circuit
-            .apply_custom(sqrt_x.clone(), &[0])
-            .apply_custom(sqrt_x.clone(), &[0]);
-        circuit
-    };
-    print_circuit(&display_circuit);
-    results.push(benchmark_circuit("√X gate (1 qubit)", builder));
-
-    let mut display = {
-        let mut circuit = QuantumCircuit::new(1);
-        circuit
-            .apply_custom(sqrt_x.clone(), &[0])
-            .apply_custom(sqrt_x.clone(), &[0]);
-        circuit
-    };
-    display.compute();
-    println!("{}", display);
-    println!("(Two √X gates should equal X, so |0⟩ becomes |1⟩)\n");
-}
-
diff --git a/tester/kernels.c b/tester/kernels.c
new file mode 100644
index 0000000..91a5cba
--- /dev/null
+++ b/tester/kernels.c
@@ -0,0 +1,53 @@
+#include "tests.h"
+
+static const double R2 = 0.7071067811865476;
+
+static void build_fusion(struct PsiQuantumCircuit* c)
+{
+	psi_apply_h(c, 0);
+	psi_apply_z(c, 0);
+	psi_apply_h(c, 0);
+	psi_apply_x(c, 1);
+	psi_apply_x(c, 1);
+}
+
+void run_kernel_tests(void)
+{
+	psi_test_section("Kernel batching and fusion");
+
+	struct PsiQuantumCircuit hh = psi_new_quantum_circuit(1);
+	psi_apply_h(&hh, 0);
+	psi_apply_h(&hh, 0);
+	struct PsiComplex hh_exp[] = { psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0) };
+	psi_check_circuit("H then H = identity", &hh, hh_exp, 2);
+	psi_free_quantum_circuit(&hh);
+
+	struct PsiKernel a =
+			psi_new_kernel("H",
+			               psi_matrix(2, 2, psi_new_complex(R2, 0.0), psi_new_complex(R2, 0.0),
+			                          psi_new_complex(R2, 0.0), psi_new_complex(-R2, 0.0)),
+			               (size_t[]){ 0 }, 1);
+	struct PsiKernel b = psi_clone_kernel(a);
+	psi_test_check(psi_kernels_can_fuse(a, b), "adjacent single-qubit kernels fuse");
+	psi_free_kernel(&a);
+	psi_free_kernel(&b);
+
+	struct PsiKernelBatch batch = psi_new_kernel_batch(1);
+	psi_add_kernel(
+			&batch,
+			psi_new_kernel("H",
+			               psi_matrix(2, 2, psi_new_complex(R2, 0.0), psi_new_complex(R2, 0.0),
+			                          psi_new_complex(R2, 0.0), psi_new_complex(-R2, 0.0)),
+			               (size_t[]){ 0 }, 1));
+	psi_add_kernel(
+			&batch,
+			psi_new_kernel("H",
+			               psi_matrix(2, 2, psi_new_complex(R2, 0.0), psi_new_complex(R2, 0.0),
+			                          psi_new_complex(R2, 0.0), psi_new_complex(-R2, 0.0)),
+			               (size_t[]){ 0 }, 1));
+	psi_optimize_kernel_batch(&batch);
+	psi_test_check(batch.count == 1, "batch fuses two H kernels into one");
+	psi_free_kernel_batch(&batch);
+
+	psi_check_runtimes_agree("Runtimes agree on fusion circuit", 2, build_fusion);
+}
diff --git a/tester/kernels.rs b/tester/kernels.rs
deleted file mode 100644
index b42b181..0000000
--- a/tester/kernels.rs
+++ /dev/null
@@ -1,420 +0,0 @@
-use crate::common::{print_section, states_equal, BenchmarkResult, CircuitCases};
-use psi::{QuantumCircuit, Runtime, RuntimeConfig};
-use std::f64::consts::PI;
-use std::time::Instant;
-
-pub fn run_all(results: &mut Vec<BenchmarkResult>) {
-    println!("═══════════════════════════════════════════════════════════════");
-    println!("                    KERNEL BATCHING TESTS");
-    println!("═══════════════════════════════════════════════════════════════\n");
-
-    test_kernel_fusion(results);
-    test_batched_vs_basic(results);
-    test_batched_large_circuits(results);
-    test_structure_aware(results);
-    test_composable_runtime(results);
-}
-
-pub fn test_kernel_fusion(results: &mut Vec<BenchmarkResult>) {
-    print_section("Kernel Fusion Test");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(2);
-        circuit.h(0).t(0).s(0).x(0).h(1).z(1);
-        circuit
-    };
-
-    let circuit = builder();
-    let batch = Runtime::build_kernel_batch(2, circuit.operations());
-    let original_count = batch.len();
-    println!("Original kernels: {}", original_count);
-    for (i, k) in batch.kernels().iter().enumerate() {
-        println!("  {}: {} on {:?}", i, k.name, k.targets);
-    }
-
-    let mut optimized_batch = Runtime::build_kernel_batch(2, circuit.operations());
-    optimized_batch.optimize();
-    let optimized_count = optimized_batch.len();
-    println!("\nOptimized kernels: {}", optimized_count);
-    for (i, k) in optimized_batch.kernels().iter().enumerate() {
-        println!("  {}: {} on {:?}", i, k.name, k.targets);
-    }
-
-    let reduction = ((original_count - optimized_count) as f64 / original_count as f64) * 100.0;
-    println!(
-        "\nKernel reduction: {} → {} ({:.0}% fewer)",
-        original_count, optimized_count, reduction
-    );
-
-    let mut basic = builder();
-    let start = Instant::now();
-    basic.compute_with(Runtime::BasicRT);
-    let basic_time = start.elapsed();
-
-    let mut batched = builder();
-    let start = Instant::now();
-    batched.compute_with(Runtime::BatchedRT);
-    let batched_time = start.elapsed();
-
-    let match_result = states_equal(basic.state(), batched.state());
-    println!("Results match: {}\n", if match_result { "✓" } else { "✗" });
-
-    results.push(BenchmarkResult {
-        name: format!("Fusion ({}→{} kernels)", original_count, optimized_count),
-        basic_time,
-        mt_time: batched_time,
-        results_match: match_result,
-    });
-
-    let fusion_heavy = || {
-        let mut circuit = QuantumCircuit::new(1);
-        circuit.h(0).t(0).s(0).x(0).y(0).z(0).h(0).t(0);
-        circuit
-    };
-
-    let circuit2 = fusion_heavy();
-    let batch2 = Runtime::build_kernel_batch(1, circuit2.operations());
-    let orig2 = batch2.len();
-    let mut opt_batch2 = Runtime::build_kernel_batch(1, circuit2.operations());
-    opt_batch2.optimize();
-    let opt2 = opt_batch2.len();
-
-    let mut basic2 = fusion_heavy();
-    let start = Instant::now();
-    basic2.compute_with(Runtime::BasicRT);
-    let basic_time2 = start.elapsed();
-
-    let mut batched2 = fusion_heavy();
-    let start = Instant::now();
-    batched2.compute_with(Runtime::BatchedRT);
-    let batched_time2 = start.elapsed();
-
-    let match2 = states_equal(basic2.state(), batched2.state());
-
-    results.push(BenchmarkResult {
-        name: format!("Heavy fusion ({}→{} kernels)", orig2, opt2),
-        basic_time: basic_time2,
-        mt_time: batched_time2,
-        results_match: match2,
-    });
-}
-
-pub fn test_batched_vs_basic(results: &mut Vec<BenchmarkResult>) {
-    print_section("Batched vs Basic Runtime Comparison");
-
-    let test_cases: CircuitCases = vec![
-        (
-            "Bell State",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(2);
-                c.h(0).cnot(0, 1);
-                c
-            }),
-        ),
-        (
-            "GHZ State",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(3);
-                c.h(0).cnot(0, 1).cnot(0, 2);
-                c
-            }),
-        ),
-        (
-            "Rotation Chain",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(3);
-                c.rx(0, PI / 4.0)
-                    .ry(0, PI / 4.0)
-                    .rz(0, PI / 4.0)
-                    .rx(1, PI / 3.0)
-                    .ry(1, PI / 3.0);
-                c
-            }),
-        ),
-        (
-            "Mixed Gates",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(4);
-                c.h(0).h(1).h(2).h(3).cnot(0, 1).cnot(2, 3).cz(1, 2);
-                c
-            }),
-        ),
-    ];
-
-    for (name, builder) in test_cases {
-        let mut basic = builder();
-        let start = Instant::now();
-        basic.compute_with(Runtime::BasicRT);
-        let basic_time = start.elapsed();
-
-        let mut batched = builder();
-        let start = Instant::now();
-        batched.compute_with(Runtime::BatchedRT);
-        let batched_time = start.elapsed();
-
-        let match_result = states_equal(basic.state(), batched.state());
-
-        println!(
-            "{}: Basic={:.2}μs, Batched={:.2}μs, Match={}",
-            name,
-            basic_time.as_secs_f64() * 1_000_000.0,
-            batched_time.as_secs_f64() * 1_000_000.0,
-            if match_result { "✓" } else { "✗" }
-        );
-
-        results.push(BenchmarkResult {
-            name: format!("Batched: {}", name),
-            basic_time,
-            mt_time: batched_time,
-            results_match: match_result,
-        });
-    }
-    println!();
-}
-
-pub fn test_batched_large_circuits(results: &mut Vec<BenchmarkResult>) {
-    print_section("Batched Runtime on Large Circuits");
-
-    let sizes = [8, 10, 12];
-
-    for &n in &sizes {
-        let builder = || {
-            let mut circuit = QuantumCircuit::new(n);
-            for i in 0..n {
-                circuit.h(i);
-            }
-            for i in 0..(n - 1) {
-                circuit.cnot(i, i + 1);
-            }
-            circuit
-        };
-
-        let mut basic_mt = builder();
-        let start = Instant::now();
-        basic_mt.compute_with(Runtime::BasicRTMT);
-        let basic_mt_time = start.elapsed();
-
-        let mut batched_mt = builder();
-        let start = Instant::now();
-        batched_mt.compute_with(Runtime::BatchedRTMT);
-        let batched_mt_time = start.elapsed();
-
-        let match_result = states_equal(basic_mt.state(), batched_mt.state());
-
-        println!(
-            "{}-qubit: BasicRTMT={:.3}ms, BatchedRTMT={:.3}ms, Match={}",
-            n,
-            basic_mt_time.as_secs_f64() * 1000.0,
-            batched_mt_time.as_secs_f64() * 1000.0,
-            if match_result { "✓" } else { "✗" }
-        );
-
-        results.push(BenchmarkResult {
-            name: format!("{}-qubit batched", n),
-            basic_time: basic_mt_time,
-            mt_time: batched_mt_time,
-            results_match: match_result,
-        });
-    }
-    println!();
-}
-
-pub fn test_structure_aware(results: &mut Vec<BenchmarkResult>) {
-    print_section("Structure-Aware Kernel Optimisation");
-
-    let commute_test = || {
-        let mut c = QuantumCircuit::new(3);
-        c.t(0).h(1).t(0).h(2).s(0).t(1).rz(0, PI / 4.0);
-        c
-    };
-
-    let circuit = commute_test();
-    let mut batch = Runtime::build_structure_aware_batch(3, circuit.operations());
-    let original = batch.len();
-    println!("Original operations: {}", original);
-    for (i, k) in batch.kernels().iter().enumerate() {
-        println!("  {}: {} on {:?} ({:?})", i, k.name, k.targets, k.gate_type);
-    }
-
-    batch.optimise();
-    let optimised = batch.len();
-    println!("\nAfter optimisation: {}", optimised);
-    for (i, k) in batch.kernels().iter().enumerate() {
-        println!("  {}: {} on {:?}", i, k.name, k.targets);
-    }
-
-    println!("\nExecution layers: {}", batch.num_layers());
-    for (i, layer) in batch.layers().iter().enumerate() {
-        let names: Vec<_> = layer.kernels.iter().map(|k| k.name.as_str()).collect();
-        println!("  Layer {}: {:?}", i, names);
-    }
-
-    let stats = batch.stats();
-    println!("\nStats: {}", stats);
-
-    let mut basic = commute_test();
-    let start = Instant::now();
-    basic.compute_with(Runtime::BasicRT);
-    let basic_time = start.elapsed();
-
-    let mut sa = commute_test();
-    let start = Instant::now();
-    sa.compute_with(Runtime::StructureAwareRT);
-    let sa_time = start.elapsed();
-
-    let match_result = states_equal(basic.state(), sa.state());
-    println!(
-        "\nBasic={:.2}μs, StructureAware={:.2}μs, Match={}",
-        basic_time.as_secs_f64() * 1_000_000.0,
-        sa_time.as_secs_f64() * 1_000_000.0,
-        if match_result { "✓" } else { "✗" }
-    );
-
-    results.push(BenchmarkResult {
-        name: format!("SA: Commuting ({}→{})", original, optimised),
-        basic_time,
-        mt_time: sa_time,
-        results_match: match_result,
-    });
-
-    println!();
-    print_section("Structure-Aware vs Other Runtimes");
-
-    let test_cases: CircuitCases = vec![
-        (
-            "Diagonal-heavy (5q)",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(5);
-                for q in 0..5 {
-                    c.t(q).s(q).rz(q, PI / 4.0).t(q);
-                }
-                c
-            }),
-        ),
-        (
-            "Interleaved (4q)",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(4);
-                c.h(0).h(1).h(2).h(3);
-                c.t(0).t(1).t(2).t(3);
-                c.cnot(0, 1).cnot(2, 3);
-                c.s(0).s(1).s(2).s(3);
-                c
-            }),
-        ),
-        (
-            "Deep rotation (3q)",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(3);
-                for _ in 0..5 {
-                    for q in 0..3 {
-                        c.rx(q, PI / 8.0).ry(q, PI / 8.0).rz(q, PI / 8.0);
-                    }
-                }
-                c
-            }),
-        ),
-    ];
-
-    for (name, builder) in test_cases {
-        let mut batched = builder();
-        let start = Instant::now();
-        batched.compute_with(Runtime::BatchedRT);
-        let batched_time = start.elapsed();
-
-        let mut sa = builder();
-        let start = Instant::now();
-        sa.compute_with(Runtime::StructureAwareRT);
-        let sa_time = start.elapsed();
-
-        let match_result = states_equal(batched.state(), sa.state());
-
-        let speedup = batched_time.as_secs_f64() / sa_time.as_secs_f64();
-        println!(
-            "{}: Batched={:.2}μs, SA={:.2}μs, Speedup={:.2}x, Match={}",
-            name,
-            batched_time.as_secs_f64() * 1_000_000.0,
-            sa_time.as_secs_f64() * 1_000_000.0,
-            speedup,
-            if match_result { "✓" } else { "✗" }
-        );
-
-        results.push(BenchmarkResult {
-            name: format!("SA: {}", name),
-            basic_time: batched_time,
-            mt_time: sa_time,
-            results_match: match_result,
-        });
-    }
-    println!();
-}
-
-pub fn test_composable_runtime(results: &mut Vec<BenchmarkResult>) {
-    print_section("Composable Runtime Configurations");
-
-    let builder = || {
-        let mut c = QuantumCircuit::new(6);
-        for q in 0..6 {
-            c.h(q).t(q).s(q);
-        }
-        for q in 0..5 {
-            c.cnot(q, q + 1);
-        }
-        for q in 0..6 {
-            c.rx(q, PI / 4.0).rz(q, PI / 4.0);
-        }
-        c
-    };
-
-    let configs: Vec<(&str, RuntimeConfig)> = vec![
-        ("Basic", RuntimeConfig::new()),
-        ("Batched", RuntimeConfig::new().batched()),
-        ("SIMD", RuntimeConfig::new().simd()),
-        ("Batched+SIMD", RuntimeConfig::new().batched().simd()),
-        ("SA+SIMD", RuntimeConfig::new().structure_aware().simd()),
-        (
-            "SA+SIMD+Parallel",
-            RuntimeConfig::new().structure_aware().simd().parallel(),
-        ),
-        ("Optimal", Runtime::optimal()),
-    ];
-
-    let mut reference = builder();
-    reference.compute_with(Runtime::BasicRT);
-    let ref_state = reference.state().clone();
-
-    println!("Testing 6-qubit circuit with different runtime configurations:\n");
-
-    for (name, config) in &configs {
-        let mut circuit = builder();
-        let start = Instant::now();
-        circuit.compute_with_config(*config);
-        let time = start.elapsed();
-
-        let match_result = states_equal(&ref_state, circuit.state());
-
-        println!(
-            "{:20} : {:.2}μs, Match={}",
-            name,
-            time.as_secs_f64() * 1_000_000.0,
-            if match_result { "✓" } else { "✗" }
-        );
-
-        results.push(BenchmarkResult {
-            name: format!("Config: {}", name),
-            basic_time: time,
-            mt_time: time,
-            results_match: match_result,
-        });
-    }
-
-    println!("\nConfiguration Display Examples:");
-    println!("  {}", RuntimeConfig::new());
-    println!("  {}", RuntimeConfig::new().batched().simd());
-    println!(
-        "  {}",
-        RuntimeConfig::new().structure_aware().simd().parallel()
-    );
-    println!("  {}", Runtime::optimal());
-    println!();
-}
diff --git a/tester/main.c b/tester/main.c
index 10d0b7d..c516568 100644
--- a/tester/main.c
+++ b/tester/main.c
@@ -1,250 +1,41 @@
 #include <stdio.h>
-#include <stdlib.h>
+#include <string.h>
 
-#include "psi.h"
+#include "tests.h"
 
-int main(void)
+static bool has_arg(int argc, char** argv, const char* name)
 {
-	printf("psi %s\n", psi_version());
+	for (int i = 1; i < argc; i++)
+		if (strcmp(argv[i], name) == 0)
+			return true;
 
-	struct PsiComplex a = psi_new_complex(1.0, 2.0);
-	struct PsiComplex b = psi_new_complex(3.0, -1.0);
-
-	struct PsiComplex sum = psi_add_complex(a, b);
-	struct PsiComplex product = psi_mul_complex(a, b);
-
-	printf("a       = %g + %gi\n", a.real, a.imaginary);
-	printf("b       = %g + %gi\n", b.real, b.imaginary);
-	printf("a + b   = %g + %gi\n", sum.real, sum.imaginary);
-	printf("a * b   = %g + %gi\n", product.real, product.imaginary);
-	printf("|a|     = %g\n", psi_abs_complex(a));
-
-	struct PsiVector u = psi_column_vector(
-		psi_new_complex(1.0, 0.0),
-		psi_new_complex(2.0, 0.0),
-		psi_new_complex(3.0, 0.0));
-	struct PsiVector v = psi_column_vector(
-		psi_new_complex(4.0, 0.0),
-		psi_new_complex(5.0, 0.0),
-		psi_new_complex(6.0, 0.0));
-
-	struct PsiComplex dot = psi_dot_vector(u, v);
-	struct PsiComplex norm = psi_norm_vector(u);
-
-	printf("u . v   = %g + %gi\n", dot.real, dot.imaginary);
-	printf("|u|     = %g + %gi\n", norm.real, norm.imaginary);
-
-	psi_free_vector(&u);
-	psi_free_vector(&v);
-
-	struct PsiMatrix m = psi_matrix(2, 2,
-		psi_new_complex(1.0, 0.0), psi_new_complex(2.0, 0.0),
-		psi_new_complex(3.0, 0.0), psi_new_complex(4.0, 0.0));
-	struct PsiMatrix mm = psi_dot_matrix(m, m);
-
-	printf("m*m     = [%g %g; %g %g]\n",
-		psi_get_matrix(mm, 0, 0).real, psi_get_matrix(mm, 0, 1).real,
-		psi_get_matrix(mm, 1, 0).real, psi_get_matrix(mm, 1, 1).real);
-
-	struct PsiVector x = psi_column_vector(
-		psi_new_complex(1.0, 0.0),
-		psi_new_complex(1.0, 0.0));
-	struct PsiVector mx = psi_mul_vector_matrix(x, m);
-
-	printf("m*x     = [%g, %g]\n", mx.data[0].real, mx.data[1].real);
-
-	psi_free_vector(&x);
-	psi_free_vector(&mx);
-	psi_free_matrix(&m);
-	psi_free_matrix(&mm);
-
-	double inv_sqrt2 = 0.7071067811865475;
-	struct PsiMatrix h = psi_matrix(2, 2,
-		psi_new_complex(inv_sqrt2, 0.0), psi_new_complex(inv_sqrt2, 0.0),
-		psi_new_complex(inv_sqrt2, 0.0), psi_new_complex(-inv_sqrt2, 0.0));
-	struct PsiQuantumGate gate = psi_new_quantum_gate_from_matrix("H", h);
-
-	struct PsiVector s0 = psi_new_state_0();
-	struct PsiVector hs0 = psi_mul_vector_matrix(s0, gate.matrix);
-
-	printf("gate    = %s (%zu qubit)\n", gate.name, gate.num_qubits);
-	printf("H|0>    = [%g, %g]\n", hs0.data[0].real, hs0.data[1].real);
-
-	psi_free_vector(&s0);
-	psi_free_vector(&hs0);
-	psi_free_quantum_gate(&gate);
-
-	const char *names[] = { "q0", "q1" };
-	struct PsiQuantumRegister reg = psi_new_quantum_register("bell", names, 2);
-
-	struct PsiQuantumGate hadamard = psi_hadamard_gate();
-	struct PsiQuantumGate cnot = psi_cnot_gate();
-
-	size_t t0[] = { 0 };
-	size_t t01[] = { 0, 1 };
-	psi_apply_gate(&reg, hadamard, t0, 1);
-	psi_apply_gate(&reg, cnot, t01, 2);
-
-	printf("bell    = [%g, %g, %g, %g]\n",
-		reg.state_vector.data[0].real, reg.state_vector.data[1].real,
-		reg.state_vector.data[2].real, reg.state_vector.data[3].real);
-
-	psi_free_quantum_gate(&hadamard);
-	psi_free_quantum_gate(&cnot);
-	psi_free_quantum_register(&reg);
-
-	struct PsiCompositeGateOp bell_ops[] = {
-		{ PSI_OP_H, { 0 }, 1 },
-		{ PSI_OP_CNOT, { 0, 1 }, 2 },
-	};
-	struct PsiCustomGate bell = psi_new_custom_gate_from_composite("BELL", 2, bell_ops, 2);
-	struct PsiQuantumGate bell_gate = psi_to_quantum_gate(bell);
-
-	struct PsiQuantumRegister reg2 = psi_new_quantum_register("bell2", names, 2);
-	size_t t01b[] = { 0, 1 };
-	psi_apply_gate(&reg2, bell_gate, t01b, 2);
-
-	printf("custom  = %s -> [%g, %g, %g, %g]\n", bell_gate.name,
-		reg2.state_vector.data[0].real, reg2.state_vector.data[1].real,
-		reg2.state_vector.data[2].real, reg2.state_vector.data[3].real);
-
-	psi_free_custom_gate(&bell);
-	psi_free_quantum_gate(&bell_gate);
-	psi_free_quantum_register(&reg2);
-
-	const char *cnames[] = { "c0", "c1" };
-	struct PsiClassicalRegister creg = psi_new_classical_register("c", cnames, 2);
-	creg.bits[1].state = true;
-
-	printf("creg    = %s=%d %s=%d\n",
-		creg.bits[0].name, creg.bits[0].state,
-		creg.bits[1].name, creg.bits[1].state);
-
-	psi_free_classical_register(&creg);
-
-	struct PsiQuantumCircuit circuit = psi_new_quantum_circuit(3);
-	psi_apply_h(&circuit, 0);
-	psi_apply_cnot(&circuit, 0, 1);
-	psi_apply_rx(&circuit, 2, 1.5707963267948966);
-	psi_apply_ccnot(&circuit, 0, 1, 2);
-
-	struct PsiCompositeGateOp cops[] = { { PSI_OP_X, { 0 }, 1 } };
-	struct PsiCustomGate xg = psi_new_custom_gate_from_composite("MYX", 1, cops, 1);
-	size_t ct[] = { 2 };
-	psi_apply_custom(&circuit, xg, ct, 1);
-
-	psi_measure_all(&circuit);
+	return false;
+}
 
-	printf("circuit = %zu qubits, %zu classical, %zu ops\n",
-		circuit.num_qubits, circuit.num_classical, circuit.operation_count);
-	for (size_t i = 0; i < circuit.operation_count; i++)
+int main(int argc, char** argv)
+{
+	if (has_arg(argc, argv, "help") || has_arg(argc, argv, "--help") || has_arg(argc, argv, "-h"))
 	{
-		size_t tc;
-		const size_t *targets = psi_gate_op_quantum_targets(&circuit.operations[i], &tc);
-		printf("  %zu: %s on", i, psi_gate_op_name(circuit.operations[i]));
-		for (size_t j = 0; j < tc; j++)
-			printf(" q%zu", targets[j]);
-		printf("\n");
+		printf("Usage: tester [clifford|non-clifford|custom|kernels|simd|noise|all]\n");
+		return 0;
 	}
 
-	psi_free_quantum_circuit(&circuit);
-
-	size_t kt[] = { 0 };
-	struct PsiKernelBatch batch = psi_new_kernel_batch(1);
-	struct PsiQuantumGate hg = psi_hadamard_gate();
-	psi_add_kernel(&batch, psi_new_kernel("H", psi_clone_matrix(hg.matrix), kt, 1));
-	psi_add_kernel(&batch, psi_new_kernel("H", psi_clone_matrix(hg.matrix), kt, 1));
-	psi_free_quantum_gate(&hg);
-
-	size_t before = batch.count;
-	psi_optimize_kernel_batch(&batch);
-
-	struct PsiVector state = psi_new_state_0();
-	psi_execute_kernel_batch(batch, &state);
-
-	printf("kernels = %zu -> %zu (fused), HH|0> = [%g, %g]\n",
-		before, batch.count, state.data[0].real, state.data[1].real);
-
-	psi_free_vector(&state);
-	psi_free_kernel_batch(&batch);
-
-	struct PsiStructureAwareBatch sa = psi_new_structure_aware_batch(2);
-	struct PsiQuantumGate hg2 = psi_hadamard_gate();
-	struct PsiQuantumGate xg2 = psi_pauli_x_gate();
-	size_t q0[] = { 0 };
-	size_t q1[] = { 1 };
-	psi_add_structure_aware_kernel(&sa, psi_new_kernel("H", psi_clone_matrix(hg2.matrix), q0, 1));
-	psi_add_structure_aware_kernel(&sa, psi_new_kernel("X", psi_clone_matrix(xg2.matrix), q1, 1));
-	psi_add_structure_aware_kernel(&sa, psi_new_kernel("H", psi_clone_matrix(hg2.matrix), q0, 1));
-	psi_free_quantum_gate(&hg2);
-	psi_free_quantum_gate(&xg2);
-
-	psi_optimize_structure_aware_batch(&sa);
-	struct PsiKernelStats stats = psi_structure_aware_batch_stats(sa);
-
-	struct PsiQuantumRegister r = psi_new_quantum_register("r", names, 2);
-	struct PsiVector sast = psi_clone_vector(r.state_vector);
-	psi_free_quantum_register(&r);
-	psi_execute_structure_aware_batch_layered(sa, &sast);
+	bool all = argc < 2 || has_arg(argc, argv, "all");
 
-	printf("sa      = %zu kernels, %zu layers -> [%g, %g, %g, %g]\n",
-		stats.total_kernels, stats.execution_layers,
-		sast.data[0].real, sast.data[1].real, sast.data[2].real, sast.data[3].real);
+	printf("psi %s test suite\n", psi_version());
 
-	psi_free_vector(&sast);
-	psi_free_structure_aware_batch(&sa);
+	if (all || has_arg(argc, argv, "clifford"))
+		run_clifford_tests();
+	if (all || has_arg(argc, argv, "non-clifford"))
+		run_non_clifford_tests();
+	if (all || has_arg(argc, argv, "custom"))
+		run_custom_tests();
+	if (all || has_arg(argc, argv, "kernels"))
+		run_kernel_tests();
+	if (all || has_arg(argc, argv, "simd"))
+		run_simd_tests();
+	if (all || has_arg(argc, argv, "noise"))
+		run_noise_tests();
 
-	struct PsiQuantumCircuit ghz = psi_new_quantum_circuit(3);
-	psi_apply_h(&ghz, 0);
-	psi_apply_cnot(&ghz, 0, 1);
-	psi_apply_cnot(&ghz, 1, 2);
-
-	const struct PsiVector *ghz_state = psi_compute_circuit_with(&ghz, PSI_RUNTIME_STRUCTURE_AWARE);
-	printf("ghz     =");
-	for (size_t i = 0; i < ghz_state->size; i++)
-		printf(" %g", ghz_state->data[i].real);
-	printf("\n");
-
-	double probs[8];
-	psi_circuit_probabilities(&ghz, probs);
-	printf("p(000)  = %g, p(111) = %g\n", probs[0], probs[7]);
-
-	char abuf[64];
-	char pbuf[64];
-	printf("amp     = |000> %s, p|111> %s\n",
-		psi_format_amplitude(ghz_state->data[0], abuf, sizeof abuf),
-		psi_format_probability(probs[7], pbuf, sizeof pbuf));
-
-	psi_free_quantum_circuit(&ghz);
-
-	struct PsiDensityMatrix dm = psi_new_density_matrix(1);
-	printf("dm      = purity %g (pure? %d)\n", psi_purity_density_matrix(dm), psi_is_pure_density_matrix(dm, 1e-10));
-
-	struct PsiNoiseChannel bf = psi_bit_flip_channel(0.25);
-	psi_apply_noise_channel(&dm, bf, 0);
-
-	double dprobs[2];
-	psi_density_matrix_probabilities(dm, dprobs);
-	printf("noisy   = p [%g, %g], purity %g\n", dprobs[0], dprobs[1], psi_purity_density_matrix(dm));
-
-	psi_free_noise_channel(&bf);
-	psi_free_density_matrix(&dm);
-
-	struct PsiQuantumCircuit viz = psi_new_quantum_circuit(3);
-	psi_apply_h(&viz, 0);
-	psi_apply_cnot(&viz, 0, 1);
-	psi_apply_rx(&viz, 2, 1.5707963267948966);
-	psi_apply_ccnot(&viz, 0, 1, 2);
-	psi_measure(&viz, 0, 0);
-
-	char *diagram = psi_render_circuit_vertical(&viz);
-	printf("\n%s\n", diagram);
-	free(diagram);
-
-	char *hdiagram = psi_render_circuit_horizontal(&viz);
-	printf("%s\n", hdiagram);
-	free(hdiagram);
-
-	psi_free_quantum_circuit(&viz);
+	return psi_test_summary();
 }
-
diff --git a/tester/main.rs b/tester/main.rs
deleted file mode 100644
index ecebd34..0000000
--- a/tester/main.rs
+++ /dev/null
@@ -1,101 +0,0 @@
-mod benchmarks;
-mod clifford;
-mod common;
-mod custom_gates;
-mod kernels;
-mod noise;
-mod non_clifford;
-mod simd;
-
-use common::{print_benchmark_table, print_summary, BenchmarkResult};
-use std::env;
-
-fn print_header() {
-    println!("═══════════════════════════════════════════════════════════════");
-    println!("                    PSI Quantum Simulator");
-    println!("═══════════════════════════════════════════════════════════════\n");
-}
-
-fn print_usage() {
-    println!("Usage: tester [OPTIONS]");
-    println!();
-    println!("Options:");
-    println!("  all          Run all tests (default)");
-    println!("  clifford     Run Clifford gate tests only");
-    println!("  non-clifford Run non-Clifford gate tests only");
-    println!("  custom       Run custom gate tests only");
-    println!("  kernels      Run kernel batching tests only");
-    println!("  simd         Run SIMD acceleration tests only");
-    println!("  noise        Run noise channel tests only");
-    println!("  bench        Run benchmark tests only");
-    println!("  help         Show this help message");
-    println!();
-    println!("Examples:");
-    println!("  tester                   # Run all tests");
-    println!("  tester clifford          # Run only Clifford gate tests");
-    println!("  tester non-clifford      # Run only rotation/parametric gate tests");
-    println!("  tester kernels           # Run only kernel batching tests");
-    println!("  tester simd              # Run only SIMD tests");
-    println!("  tester noise             # Run only noise channel tests");
-    println!("  tester custom bench      # Run custom gates and benchmarks");
-}
-
-fn main() {
-    let args: Vec<String> = env::args().skip(1).collect();
-
-    if args
-        .iter()
-        .any(|a| a == "help" || a == "--help" || a == "-h")
-    {
-        print_usage();
-        return;
-    }
-
-    print_header();
-
-    let mut results: Vec<BenchmarkResult> = Vec::new();
-
-    let run_all = args.is_empty() || args.iter().any(|a| a == "all");
-    let run_clifford = run_all || args.iter().any(|a| a == "clifford");
-    let run_non_clifford = run_all || args.iter().any(|a| a == "non-clifford");
-    let run_custom = run_all || args.iter().any(|a| a == "custom");
-    let run_kernels = run_all || args.iter().any(|a| a == "kernels");
-    let run_simd = run_all || args.iter().any(|a| a == "simd");
-    let run_noise = run_all || args.iter().any(|a| a == "noise");
-    let run_bench = run_all || args.iter().any(|a| a == "bench");
-
-    if run_clifford {
-        clifford::run_all(&mut results);
-    }
-    
-
-    if run_non_clifford {
-        non_clifford::run_all(&mut results);
-    }
-
-    if run_custom {
-        custom_gates::run_all(&mut results);
-    }
-
-    if run_kernels {
-        kernels::run_all(&mut results);
-    }
-
-    if run_simd {
-        simd::run_all(&mut results);
-    }
-
-    if run_noise {
-        noise::run_all(&mut results);
-    }
-
-    if run_bench {
-        benchmarks::run_all(&mut results);
-    }
-    
-
-    if !results.is_empty() {
-        print_benchmark_table(&results);
-        print_summary(&results);
-    }
-}
diff --git a/tester/noise.c b/tester/noise.c
new file mode 100644
index 0000000..df791ae
--- /dev/null
+++ b/tester/noise.c
@@ -0,0 +1,49 @@
+#include "tests.h"
+
+#include <math.h>
+
+void run_noise_tests(void)
+{
+	psi_test_section("Noise channels");
+
+	struct PsiDensityMatrix pure = psi_new_density_matrix(1);
+	psi_test_check(psi_is_pure_density_matrix(pure, 1e-10), "fresh density matrix is pure");
+	psi_test_check(fabs(psi_trace_density_matrix(pure).real - 1.0) < 1e-10, "trace = 1");
+	psi_free_density_matrix(&pure);
+
+	struct PsiDensityMatrix bf = psi_new_density_matrix(1);
+	struct PsiNoiseChannel flip = psi_bit_flip_channel(0.25);
+	psi_apply_noise_channel(&bf, flip, 0);
+	double probs[2];
+	psi_density_matrix_probabilities(bf, probs);
+	psi_test_check(fabs(probs[0] - 0.75) < 1e-10 && fabs(probs[1] - 0.25) < 1e-10,
+	               "bit_flip(0.25) on |0>: p = [0.75, 0.25]");
+	psi_test_check(fabs(psi_purity_density_matrix(bf) - 0.625) < 1e-10,
+	               "bit_flip(0.25) purity = 0.625");
+	psi_free_noise_channel(&flip);
+	psi_free_density_matrix(&bf);
+
+	struct PsiComplex one[2] = { psi_new_complex(0.0, 0.0), psi_new_complex(1.0, 0.0) };
+	struct PsiDensityMatrix damp = psi_new_density_matrix_from_state(one, 2);
+	struct PsiNoiseChannel ad = psi_amplitude_damping_channel(1.0);
+	psi_apply_noise_channel(&damp, ad, 0);
+	double dprobs[2];
+	psi_density_matrix_probabilities(damp, dprobs);
+	psi_test_check(fabs(dprobs[0] - 1.0) < 1e-10, "amplitude_damping(1.0) fully decays |1> -> |0>");
+	psi_free_noise_channel(&ad);
+	psi_free_density_matrix(&damp);
+
+	struct PsiComplex zero[2] = { psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0) };
+	struct PsiDensityMatrix fid = psi_new_density_matrix_from_state(zero, 2);
+	psi_test_check(fabs(psi_fidelity_density_matrix(fid, zero) - 1.0) < 1e-10,
+	               "fidelity(|0><0|, |0>) = 1");
+	psi_free_density_matrix(&fid);
+
+	struct PsiDensityMatrix depo = psi_new_density_matrix(1);
+	struct PsiNoiseChannel dc = psi_depolarising_channel(0.3);
+	psi_apply_noise_channel(&depo, dc, 0);
+	psi_test_check(fabs(psi_trace_density_matrix(depo).real - 1.0) < 1e-10,
+	               "depolarising preserves trace");
+	psi_free_noise_channel(&dc);
+	psi_free_density_matrix(&depo);
+}
diff --git a/tester/noise.rs b/tester/noise.rs
deleted file mode 100644
index 02a4ba3..0000000
--- a/tester/noise.rs
+++ /dev/null
@@ -1,172 +0,0 @@
-use crate::common::{print_section, BenchmarkResult};
-use psi::{
-    complex, DensityMatrix, NoiseChannel, QuantumCircuit, Runtime, Vector,
-};
-use std::time::Instant;
-
-pub fn run_all(results: &mut Vec<BenchmarkResult>) {
-    println!("═══════════════════════════════════════════════════════════════");
-    println!("                    NOISE CHANNEL TESTS");
-    println!("═══════════════════════════════════════════════════════════════\n");
-
-    test_density_matrix_basics(results);
-    test_noise_channels(results);
-    test_noisy_circuit(results);
-}
-
-pub fn test_density_matrix_basics(results: &mut Vec<BenchmarkResult>) {
-    print_section("Density Matrix Basics");
-
-    let dm = DensityMatrix::new(2);
-    println!("Initial |00⟩ state:");
-    println!("{}", dm);
-
-    let mut circuit = QuantumCircuit::new(2);
-    circuit.h(0).cnot(0, 1);
-    circuit.compute_with(Runtime::BasicRT);
-    let state = circuit.state();
-
-    let state_vec: Vec<_> = (0..state.size())
-        .map(|i| state.get(i))
-        .collect();
-
-    let dm_bell = DensityMatrix::from_state_vector(&state_vec);
-    println!("Bell state |Φ+⟩:");
-    println!("{}", dm_bell);
-    println!("Full matrix:");
-    println!("{:?}", dm_bell);
-
-    let is_pure = dm_bell.is_pure(1e-10);
-    println!("Purity check: {}\n", if is_pure { "✓ Pure" } else { "✗ Mixed" });
-
-    results.push(BenchmarkResult {
-        name: "DM: Bell state".to_string(),
-        basic_time: std::time::Duration::from_micros(0),
-        mt_time: std::time::Duration::from_micros(0),
-        results_match: is_pure,
-    });
-}
-
-pub fn test_noise_channels(results: &mut Vec<BenchmarkResult>) {
-    print_section("Noise Channel Effects");
-
-    let channels: Vec<(&str, NoiseChannel)> = vec![
-        ("Depolarising (p=0.1)", NoiseChannel::depolarising(0.1)),
-        ("Amplitude Damping (γ=0.2)", NoiseChannel::amplitude_damping(0.2)),
-        ("Phase Damping (γ=0.2)", NoiseChannel::phase_damping(0.2)),
-        ("Bit Flip (p=0.1)", NoiseChannel::bit_flip(0.1)),
-        ("Phase Flip (p=0.1)", NoiseChannel::phase_flip(0.1)),
-        ("Bit-Phase Flip (p=0.1)", NoiseChannel::bit_phase_flip(0.1)),
-    ];
-
-    let plus_state = vec![
-        complex!(1.0 / 2.0_f64.sqrt(), 0.0),
-        complex!(1.0 / 2.0_f64.sqrt(), 0.0),
-    ];
-
-    println!("Starting with |+⟩ state: (|0⟩ + |1⟩)/√2\n");
-
-    for (name, channel) in channels {
-        let mut dm = DensityMatrix::from_state_vector(&plus_state);
-        let initial_purity = dm.purity();
-
-        let start = Instant::now();
-        dm.apply_noise_channel(&channel, 0);
-        let elapsed = start.elapsed();
-
-        let final_purity = dm.purity();
-        let fidelity = dm.fidelity_with_pure_state(&plus_state);
-
-        println!("{:30}", name);
-        println!("  Purity: {:.4} → {:.4}", initial_purity, final_purity);
-        println!("  Fidelity with |+⟩: {:.4}", fidelity);
-        println!("  Probabilities: {:?}", dm.probabilities());
-        println!("  Time: {:.2}μs\n", elapsed.as_secs_f64() * 1_000_000.0);
-
-        let purity_decreased = final_purity <= initial_purity + 1e-10;
-
-        results.push(BenchmarkResult {
-            name: format!("Noise: {}", name),
-            basic_time: elapsed,
-            mt_time: elapsed,
-            results_match: purity_decreased,
-        });
-    }
-}
-
-pub fn test_noisy_circuit(results: &mut Vec<BenchmarkResult>) {
-    print_section("Noisy Circuit Simulation");
-
-    let mut circuit = QuantumCircuit::new(2);
-    circuit.h(0).cnot(0, 1);
-    circuit.compute_with(Runtime::BasicRT);
-    let state = circuit.state();
-    let state_vec: Vec<_> = (0..state.size()).map(|i| state.get(i)).collect();
-
-    let mut dm = DensityMatrix::from_state_vector(&state_vec);
-    println!("Bell state before noise:");
-    println!("{}", dm);
-
-    let depol = NoiseChannel::depolarising(0.05);
-
-    let start = Instant::now();
-    dm.apply_noise_channel(&depol, 0);
-    dm.apply_noise_channel(&depol, 1);
-    let elapsed = start.elapsed();
-
-    println!("Bell state after 5% depolarising on both qubits:");
-    println!("{}", dm);
-
-    let fidelity = dm.fidelity_with_pure_state(&state_vec);
-    println!("Fidelity with ideal Bell state: {:.4}", fidelity);
-    println!("Time: {:.2}μs\n", elapsed.as_secs_f64() * 1_000_000.0);
-
-    let mut dm2 = DensityMatrix::from_state_vector(&state_vec);
-    let amp_damp = NoiseChannel::amplitude_damping(0.1);
-
-    dm2.apply_noise_channel(&amp_damp, 0);
-    dm2.apply_noise_channel(&amp_damp, 1);
-
-    println!("Bell state after 10% amplitude damping on both qubits:");
-    println!("{}", dm2);
-    println!("Probabilities show decay towards |00⟩: {:?}", dm2.probabilities());
-
-    results.push(BenchmarkResult {
-        name: "Noisy Bell circuit".to_string(),
-        basic_time: elapsed,
-        mt_time: elapsed,
-        results_match: fidelity > 0.8 && fidelity < 1.0,
-    });
-
-    println!();
-    print_section("T1/T2 Relaxation Simulation");
-
-    let one_state = vec![complex!(0.0, 0.0), complex!(1.0, 0.0)];
-    let mut dm_t1 = DensityMatrix::from_state_vector(&one_state);
-
-    println!("Simulating T1 decay of |1⟩ state:");
-    println!("  Initial: P(0)={:.4}, P(1)={:.4}", dm_t1.probabilities()[0], dm_t1.probabilities()[1]);
-
-    let t1_channel = NoiseChannel::amplitude_damping(0.3);
-    for step in 1..=5 {
-        dm_t1.apply_noise_channel(&t1_channel, 0);
-        println!(
-            "  Step {}: P(0)={:.4}, P(1)={:.4}, Purity={:.4}",
-            step,
-            dm_t1.probabilities()[0],
-            dm_t1.probabilities()[1],
-            dm_t1.purity()
-        );
-    }
-
-    let decayed = dm_t1.probabilities()[0] > 0.8;
-    println!("  Decay complete: {}\n", if decayed { "✓" } else { "✗" });
-
-    results.push(BenchmarkResult {
-        name: "T1 decay simulation".to_string(),
-        basic_time: std::time::Duration::from_micros(0),
-        mt_time: std::time::Duration::from_micros(0),
-        results_match: decayed,
-    });
-}
-
diff --git a/tester/non_clifford.c b/tester/non_clifford.c
new file mode 100644
index 0000000..c361fb5
--- /dev/null
+++ b/tester/non_clifford.c
@@ -0,0 +1,42 @@
+#include "tests.h"
+
+static const double R2 = 0.7071067811865476;
+static const double PI = 3.141592653589793;
+
+void run_non_clifford_tests(void)
+{
+	psi_test_section("Non-Clifford gates");
+
+	struct PsiQuantumCircuit t = psi_new_quantum_circuit(1);
+	psi_apply_x(&t, 0);
+	psi_apply_t(&t, 0);
+	struct PsiComplex t_exp[] = { psi_new_complex(0.0, 0.0), psi_new_complex(R2, R2) };
+	psi_check_circuit("T|1> = e^{i pi/4}|1>", &t, t_exp, 2);
+	psi_free_quantum_circuit(&t);
+
+	struct PsiQuantumCircuit rx = psi_new_quantum_circuit(1);
+	psi_apply_rx(&rx, 0, PI);
+	struct PsiComplex rx_exp[] = { psi_new_complex(0.0, 0.0), psi_new_complex(0.0, -1.0) };
+	psi_check_circuit("Rx(pi)|0> = -i|1>", &rx, rx_exp, 2);
+	psi_free_quantum_circuit(&rx);
+
+	struct PsiQuantumCircuit ry = psi_new_quantum_circuit(1);
+	psi_apply_ry(&ry, 0, PI);
+	struct PsiComplex ry_exp[] = { psi_new_complex(0.0, 0.0), psi_new_complex(1.0, 0.0) };
+	psi_check_circuit("Ry(pi)|0> = |1>", &ry, ry_exp, 2);
+	psi_free_quantum_circuit(&ry);
+
+	struct PsiQuantumCircuit p = psi_new_quantum_circuit(1);
+	psi_apply_x(&p, 0);
+	psi_apply_p(&p, 0, PI / 2.0);
+	struct PsiComplex p_exp[] = { psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 1.0) };
+	psi_check_circuit("P(pi/2)|1> = i|1>", &p, p_exp, 2);
+	psi_free_quantum_circuit(&p);
+
+	struct PsiQuantumCircuit sx = psi_new_quantum_circuit(1);
+	psi_apply_sx(&sx, 0);
+	psi_apply_sx(&sx, 0);
+	struct PsiComplex sx_exp[] = { psi_new_complex(0.0, 0.0), psi_new_complex(1.0, 0.0) };
+	psi_check_circuit("sqrt(X) applied twice = X", &sx, sx_exp, 2);
+	psi_free_quantum_circuit(&sx);
+}
diff --git a/tester/non_clifford.rs b/tester/non_clifford.rs
deleted file mode 100644
index 2146fb2..0000000
--- a/tester/non_clifford.rs
+++ /dev/null
@@ -1,137 +0,0 @@
-use crate::common::{benchmark_circuit, print_circuit, print_section, BenchmarkResult};
-use psi::QuantumCircuit;
-use std::f64::consts::PI;
-
-pub fn run_all(results: &mut Vec<BenchmarkResult>) {
-    println!("═══════════════════════════════════════════════════════════════");
-    println!("                   NON-CLIFFORD GATES TESTS");
-    println!("═══════════════════════════════════════════════════════════════\n");
-
-    test_fixed_gates(results);
-    test_rotation_gates(results);
-    test_phase_gates(results);
-    test_general_unitaries(results);
-    test_controlled_rotations(results);
-    test_variational_circuit(results);
-}
-
-pub fn test_fixed_gates(results: &mut Vec<BenchmarkResult>) {
-    print_section("Non-Clifford Gates: T, T†, √X, S†");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(2);
-        circuit.h(0).t(0).tdg(0).sx(1).sxdg(1).h(0).s(0).sdg(0);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("Non-Clifford fixed gates", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_rotation_gates(results: &mut Vec<BenchmarkResult>) {
-    print_section("Rotation Gates: Rx, Ry, Rz");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(3);
-        circuit
-            .rx(0, PI / 4.0)
-            .ry(1, PI / 2.0)
-            .rz(2, PI)
-            .rx(0, -PI / 4.0);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("Rotation gates (3 qubits)", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_phase_gates(results: &mut Vec<BenchmarkResult>) {
-    print_section("Phase Gate: P(θ)");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(2);
-        circuit.h(0).p(0, PI / 4.0).h(1).p(1, PI / 2.0);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("Phase gates (2 qubits)", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_general_unitaries(results: &mut Vec<BenchmarkResult>) {
-    print_section("General Unitaries: U1, U2, U3");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(3);
-        circuit
-            .u1(0, PI / 4.0)
-            .u2(1, 0.0, PI)
-            .u3(2, PI / 2.0, 0.0, PI);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("General unitaries (3 qubits)", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_controlled_rotations(results: &mut Vec<BenchmarkResult>) {
-    print_section("Controlled Rotation Gates: CRx, CRy, CRz, CP");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(4);
-        circuit
-            .x(0)
-            .crx(0, 1, PI / 2.0)
-            .x(2)
-            .cry(2, 3, PI / 4.0)
-            .crz(0, 2, PI)
-            .cp(1, 3, PI / 2.0);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit(
-        "Controlled rotations (4 qubits)",
-        builder,
-    ));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
-
-pub fn test_variational_circuit(results: &mut Vec<BenchmarkResult>) {
-    print_section("Variational Circuit (VQE-like)");
-
-    let builder = || {
-        let mut circuit = QuantumCircuit::new(3);
-        circuit.ry(0, 0.5).ry(1, 0.3).ry(2, 0.7);
-        circuit.cnot(0, 1).cnot(1, 2);
-        circuit.rx(0, 0.2).rx(1, 0.4).rx(2, 0.6);
-        circuit.cz(0, 2);
-        circuit
-    };
-
-    print_circuit(&builder());
-    results.push(benchmark_circuit("Variational circuit (3 qubits)", builder));
-
-    let mut display = builder();
-    display.compute();
-    println!("{}\n", display);
-}
diff --git a/tester/simd.c b/tester/simd.c
new file mode 100644
index 0000000..bea2f06
--- /dev/null
+++ b/tester/simd.c
@@ -0,0 +1,63 @@
+#include <stdio.h>
+
+#include "tests.h"
+
+static const double R2 = 0.7071067811865476;
+
+void run_simd_tests(void)
+{
+	psi_test_section("SIMD single-qubit kernels");
+
+	printf("  detected: %s\n", psi_simd_name(psi_detect_simd()));
+
+	struct PsiComplex h[2][2] = {
+		{ psi_new_complex(R2, 0.0), psi_new_complex(R2, 0.0) },
+		{ psi_new_complex(R2, 0.0), psi_new_complex(-R2, 0.0) },
+	};
+
+	struct PsiComplex two[4] = {
+		psi_new_complex(1.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+	};
+	psi_apply_single_qubit_gate_simd(two, h, 0, 2);
+	struct PsiComplex two_exp[] = {
+		psi_new_complex(R2, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(R2, 0.0),
+		psi_new_complex(0.0, 0.0),
+	};
+	psi_test_check(psi_amps_match(two, two_exp, 4), "SIMD H on q0 of |00>");
+
+	struct PsiComplex three[8];
+	for (size_t i = 0; i < 8; i++)
+		three[i] = psi_new_complex(0.0, 0.0);
+	three[0] = psi_new_complex(1.0, 0.0);
+	psi_apply_single_qubit_gate_simd(three, h, 0, 3);
+	struct PsiComplex three_exp[8];
+	for (size_t i = 0; i < 8; i++)
+		three_exp[i] = psi_new_complex(0.0, 0.0);
+	three_exp[0] = psi_new_complex(R2, 0.0);
+	three_exp[4] = psi_new_complex(R2, 0.0);
+	psi_test_check(psi_amps_match(three, three_exp, 8), "SIMD H on q0 of |000> (chunked)");
+
+	struct PsiComplex x[2][2] = {
+		{ psi_new_complex(0.0, 0.0), psi_new_complex(1.0, 0.0) },
+		{ psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0) },
+	};
+	struct PsiComplex flip[4] = {
+		psi_new_complex(1.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+	};
+	psi_apply_single_qubit_gate_simd(flip, x, 1, 2);
+	struct PsiComplex flip_exp[] = {
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(1.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+		psi_new_complex(0.0, 0.0),
+	};
+	psi_test_check(psi_amps_match(flip, flip_exp, 4), "SIMD X on q1 of |00> = |01>");
+}
diff --git a/tester/simd.rs b/tester/simd.rs
deleted file mode 100644
index e91dc6e..0000000
--- a/tester/simd.rs
+++ /dev/null
@@ -1,210 +0,0 @@
-use crate::common::{print_section, states_equal, BenchmarkResult, CircuitCases};
-use psi::{get_simd_info, QuantumCircuit, Runtime};
-use std::f64::consts::PI;
-use std::time::Instant;
-
-pub fn run_all(results: &mut Vec<BenchmarkResult>) {
-    println!("═══════════════════════════════════════════════════════════════");
-    println!("                    SIMD ACCELERATION TESTS");
-    println!("═══════════════════════════════════════════════════════════════\n");
-
-    println!("Detected: {}\n", get_simd_info());
-
-    test_simd_correctness(results);
-    test_simd_vs_batched(results);
-    test_simd_large_circuits(results);
-}
-
-pub fn test_simd_correctness(results: &mut Vec<BenchmarkResult>) {
-    print_section("SIMD Correctness Verification");
-
-    let test_cases: CircuitCases = vec![
-        (
-            "Bell State",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(2);
-                c.h(0).cnot(0, 1);
-                c
-            }),
-        ),
-        (
-            "GHZ-3",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(3);
-                c.h(0).cnot(0, 1).cnot(0, 2);
-                c
-            }),
-        ),
-        (
-            "Rotation Chain",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(3);
-                c.rx(0, PI / 4.0)
-                    .ry(0, PI / 4.0)
-                    .rz(0, PI / 4.0)
-                    .rx(1, PI / 3.0)
-                    .ry(1, PI / 3.0);
-                c
-            }),
-        ),
-        (
-            "Mixed Single-Qubit",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(4);
-                c.h(0).t(0).s(0).x(0).h(1).y(1).z(1).h(2).t(2).h(3).s(3);
-                c
-            }),
-        ),
-    ];
-
-    for (name, builder) in test_cases {
-        let mut basic = builder();
-        basic.compute_with(Runtime::BasicRT);
-
-        let mut simd = builder();
-        simd.compute_with(Runtime::SimdRT);
-
-        let match_result = states_equal(basic.state(), simd.state());
-
-        println!(
-            "{}: {}",
-            name,
-            if match_result {
-                "✓ Match"
-            } else {
-                "✗ MISMATCH"
-            }
-        );
-
-        results.push(BenchmarkResult {
-            name: format!("SIMD verify: {}", name),
-            basic_time: std::time::Duration::from_micros(0),
-            mt_time: std::time::Duration::from_micros(0),
-            results_match: match_result,
-        });
-    }
-    println!();
-}
-
-pub fn test_simd_vs_batched(results: &mut Vec<BenchmarkResult>) {
-    print_section("SIMD vs Batched Runtime Comparison");
-
-    let test_cases: CircuitCases = vec![
-        (
-            "Single-Qubit Heavy (6q)",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(6);
-                for q in 0..6 {
-                    c.h(q).t(q).s(q).x(q).y(q).z(q);
-                }
-                c
-            }),
-        ),
-        (
-            "Rotation Circuit (5q)",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(5);
-                for q in 0..5 {
-                    c.rx(q, PI / 4.0).ry(q, PI / 3.0).rz(q, PI / 6.0);
-                }
-                c
-            }),
-        ),
-        (
-            "Deep Single-Qubit (4q)",
-            Box::new(|| {
-                let mut c = QuantumCircuit::new(4);
-                for _ in 0..10 {
-                    for q in 0..4 {
-                        c.h(q).t(q);
-                    }
-                }
-                c
-            }),
-        ),
-    ];
-
-    for (name, builder) in test_cases {
-        let mut batched = builder();
-        let start = Instant::now();
-        batched.compute_with(Runtime::BatchedRT);
-        let batched_time = start.elapsed();
-
-        let mut simd = builder();
-        let start = Instant::now();
-        simd.compute_with(Runtime::SimdRT);
-        let simd_time = start.elapsed();
-
-        let match_result = states_equal(batched.state(), simd.state());
-
-        let speedup = batched_time.as_secs_f64() / simd_time.as_secs_f64();
-        println!(
-            "{}: Batched={:.2}μs, SIMD={:.2}μs, Speedup={:.2}x, Match={}",
-            name,
-            batched_time.as_secs_f64() * 1_000_000.0,
-            simd_time.as_secs_f64() * 1_000_000.0,
-            speedup,
-            if match_result { "✓" } else { "✗" }
-        );
-
-        results.push(BenchmarkResult {
-            name: format!("SIMD: {}", name),
-            basic_time: batched_time,
-            mt_time: simd_time,
-            results_match: match_result,
-        });
-    }
-    println!();
-}
-
-pub fn test_simd_large_circuits(results: &mut Vec<BenchmarkResult>) {
-    print_section("SIMD on Large Circuits (Multi-threaded)");
-
-    let sizes = [8, 10, 12];
-
-    for &n in &sizes {
-        let builder = || {
-            let mut circuit = QuantumCircuit::new(n);
-            for i in 0..n {
-                circuit.h(i);
-            }
-            for i in 0..(n - 1) {
-                circuit.cnot(i, i + 1);
-            }
-            for i in 0..n {
-                circuit.t(i).s(i);
-            }
-            circuit
-        };
-
-        let mut batched_mt = builder();
-        let start = Instant::now();
-        batched_mt.compute_with(Runtime::BatchedRTMT);
-        let batched_time = start.elapsed();
-
-        let mut simd_mt = builder();
-        let start = Instant::now();
-        simd_mt.compute_with(Runtime::SimdRTMT);
-        let simd_time = start.elapsed();
-
-        let match_result = states_equal(batched_mt.state(), simd_mt.state());
-
-        let speedup = batched_time.as_secs_f64() / simd_time.as_secs_f64();
-        println!(
-            "{}-qubit: BatchedMT={:.3}ms, SIMD_MT={:.3}ms, Speedup={:.2}x, Match={}",
-            n,
-            batched_time.as_secs_f64() * 1000.0,
-            simd_time.as_secs_f64() * 1000.0,
-            speedup,
-            if match_result { "✓" } else { "✗" }
-        );
-
-        results.push(BenchmarkResult {
-            name: format!("{}-qubit SIMD", n),
-            basic_time: batched_time,
-            mt_time: simd_time,
-            results_match: match_result,
-        });
-    }
-    println!();
-}
diff --git a/tester/test.c b/tester/test.c
new file mode 100644
index 0000000..9e11d8c
--- /dev/null
+++ b/tester/test.c
@@ -0,0 +1,81 @@
+#include "tests.h"
+
+#include <math.h>
+#include <stdio.h>
+
+static int tests_run = 0;
+static int tests_passed = 0;
+
+void psi_test_section(const char* name)
+{
+	printf("\n── %s ──\n", name);
+}
+
+void psi_test_check(bool ok, const char* name)
+{
+	tests_run++;
+	if (ok)
+		tests_passed++;
+
+	printf("  [%s] %s\n", ok ? "PASS" : "FAIL", name);
+}
+
+bool psi_amps_match(const struct PsiComplex* actual, const struct PsiComplex* expected, size_t n)
+{
+	for (size_t i = 0; i < n; i++)
+	{
+		if (fabs(actual[i].real - expected[i].real) > 1e-9)
+			return false;
+		if (fabs(actual[i].imaginary - expected[i].imaginary) > 1e-9)
+			return false;
+	}
+
+	return true;
+}
+
+void psi_check_circuit(const char* name, struct PsiQuantumCircuit* circuit,
+                       const struct PsiComplex* expected, size_t n)
+{
+	const struct PsiVector* state = psi_compute_circuit(circuit);
+	bool ok = state->size == n && psi_amps_match(state->data, expected, n);
+	psi_test_check(ok, name);
+}
+
+void psi_check_runtimes_agree(const char* name, size_t num_qubits,
+                              void (*build)(struct PsiQuantumCircuit*))
+{
+	enum PsiRuntime runtimes[] = {
+		PSI_RUNTIME_BASIC,
+		PSI_RUNTIME_BATCHED,
+		PSI_RUNTIME_SIMD,
+		PSI_RUNTIME_STRUCTURE_AWARE,
+	};
+
+	struct PsiQuantumCircuit base = psi_new_quantum_circuit(num_qubits);
+	build(&base);
+	const struct PsiVector* base_state = psi_compute_circuit_with(&base, PSI_RUNTIME_BASIC);
+	struct PsiVector reference = psi_clone_vector(*base_state);
+
+	bool ok = true;
+	for (size_t i = 1; i < sizeof runtimes / sizeof runtimes[0]; i++)
+	{
+		struct PsiQuantumCircuit circuit = psi_new_quantum_circuit(num_qubits);
+		build(&circuit);
+		const struct PsiVector* state = psi_compute_circuit_with(&circuit, runtimes[i]);
+		if (state->size != reference.size ||
+		    !psi_amps_match(state->data, reference.data, reference.size))
+			ok = false;
+
+		psi_free_quantum_circuit(&circuit);
+	}
+
+	psi_free_vector(&reference);
+	psi_free_quantum_circuit(&base);
+	psi_test_check(ok, name);
+}
+
+int psi_test_summary(void)
+{
+	printf("\n%d/%d checks passed\n", tests_passed, tests_run);
+	return tests_passed == tests_run ? 0 : 1;
+}
diff --git a/tester/tests.h b/tester/tests.h
new file mode 100644
index 0000000..d5262aa
--- /dev/null
+++ b/tester/tests.h
@@ -0,0 +1,22 @@
+#pragma once
+
+#include <stdbool.h>
+#include <stddef.h>
+
+#include "psi.h"
+
+void psi_test_section(const char* name);
+void psi_test_check(bool ok, const char* name);
+bool psi_amps_match(const struct PsiComplex* actual, const struct PsiComplex* expected, size_t n);
+void psi_check_circuit(const char* name, struct PsiQuantumCircuit* circuit,
+                       const struct PsiComplex* expected, size_t n);
+void psi_check_runtimes_agree(const char* name, size_t num_qubits,
+                              void (*build)(struct PsiQuantumCircuit*));
+int psi_test_summary(void);
+
+void run_clifford_tests(void);
+void run_non_clifford_tests(void);
+void run_kernel_tests(void);
+void run_simd_tests(void);
+void run_noise_tests(void);
+void run_custom_tests(void);