feat: simd + testing + formatting

53 files changed, 1902 insertions, 2425 deletions

diff --git a/.clang-format b/.clang-format
new file mode 100644
index 0000000..530613b
--- /dev/null
+++ b/.clang-format
@@ -0,0 +1,20 @@
+Language: Cpp
+BasedOnStyle: LLVM
+IndentWidth: 4
+TabWidth: 4
+UseTab: AlignWithSpaces
+ContinuationIndentWidth: 8
+ColumnLimit: 100
+BreakBeforeBraces: Allman
+PointerAlignment: Left
+AllowShortFunctionsOnASingleLine: None
+AllowShortIfStatementsOnASingleLine: Never
+AllowShortLoopsOnASingleLine: false
+AllowShortBlocksOnASingleLine: Never
+AllowShortCaseLabelsOnASingleLine: true
+IndentCaseLabels: true
+SpaceBeforeParens: ControlStatements
+SortIncludes: false
+Cpp11BracedListStyle: false
+AlignOperands: DontAlign
+AlignAfterOpenBracket: Align
diff --git a/.clang-tidy b/.clang-tidy
new file mode 100644
index 0000000..6d55a74
--- /dev/null
+++ b/.clang-tidy
@@ -0,0 +1,12 @@
+Checks: '-*,readability-identifier-naming'
+WarningsAsErrors: 'readability-identifier-naming'
+CheckOptions:
+  readability-identifier-naming.StructCase: CamelCase
+  readability-identifier-naming.StructPrefix: Psi
+  readability-identifier-naming.EnumCase: CamelCase
+  readability-identifier-naming.EnumPrefix: Psi
+  readability-identifier-naming.EnumConstantCase: UPPER_CASE
+  readability-identifier-naming.FunctionCase: lower_case
+  readability-identifier-naming.ParameterCase: lower_case
+  readability-identifier-naming.MacroDefinitionCase: UPPER_CASE
+  readability-identifier-naming.MacroDefinitionIgnoredRegexp: '^psi_.*'
diff --git a/.gitignore b/.gitignore
index f3f4346..e1cd3b6 100644
--- a/.gitignore
+++ b/.gitignore
@@ -21,9 +21,4 @@ Makefile
 *.vcxproj
 *.vcxproj.*
 
-# rust (reference sources during the port)
-debug/
-target/
-Cargo.lock
-**/*.rs.bk
 *.pdb
diff --git a/STYLE.md b/STYLE.md
new file mode 100644
index 0000000..856287f
--- /dev/null
+++ b/STYLE.md
@@ -0,0 +1,274 @@
+# psi C style guide
+
+This document describes the conventions for the C implementation of `psi`. Formatting
+is enforced by [`.clang-format`](.clang-format) and naming by [`.clang-tidy`](.clang-tidy);
+everything a tool cannot check is described here and is expected in review.
+
+The guiding idea: **idiomatic, data-driven C that reads like the Rust original in naming**,
+with Allman bracing as the one deliberate departure from typical C style.
+
+---
+
+## 1. Language & build
+
+- **C17**, compiled with `-Wall -Wextra -Werror -pedantic`. Warnings are errors; keep every
+  build clean.
+- Build system is **premake5** (`premake5 gmake && make config=release`). Public headers live
+  in `include/`, implementation in `src/`, the test runner in `tester/`.
+- Portable across x86_64 and ARM64 (the only supported architectures). Platform-specific code
+  (SIMD intrinsics) is `#ifdef`-gated with a scalar fallback.
+
+---
+
+## 2. Formatting (enforced by clang-format)
+
+Run before committing:
+
+```bash
+find include src tester -name '*.c' -o -name '*.h' | xargs clang-format -i
+```
+
+CI check (no diffs allowed):
+
+```bash
+find include src tester -name '*.c' -o -name '*.h' | xargs clang-format --dry-run -Werror
+```
+
+### Braces — Allman, everywhere
+
+Opening brace on its own line for functions, structs, enums, unions, and control blocks.
+
+```c
+struct PsiComplex
+{
+	double real;
+	double imaginary;
+};
+
+double psi_abs_complex(struct PsiComplex z)
+{
+	return sqrt(psi_norm2_complex(z));
+}
+```
+
+**Exception (tool limitation):** compound literals keep the brace on the same line as the
+type, because clang-format cannot break it:
+
+```c
+return (struct PsiComplex){
+	a.real + b.real,
+	a.imaginary + b.imaginary,
+};
+```
+
+### Single statements omit braces
+
+An `if`/`for`/`while` with a single-statement body has no braces; the statement goes on the
+next line.
+
+```c
+if (index >= row->width)
+	return;
+
+for (size_t i = 0; i < count; i++)
+	sum = psi_add_complex(sum, values[i]);
+```
+
+Use braces once a body has more than one statement.
+
+### Indentation — tabs
+
+Tabs for indentation (width 4). Alignment (e.g. wrapped arguments under an open paren) uses
+spaces, so alignment survives any tab width.
+
+### Pointers bind to the type
+
+```c
+struct PsiComplex* data;
+void psi_free_vector(struct PsiVector* v);
+const size_t* targets;
+```
+
+Not `struct PsiComplex *data`.
+
+### Line length & wrapping
+
+- 100-column limit. Long lines wrap; continuation indents two tabs, wrapped call/parameter
+  lists align under the opening paren.
+
+```c
+static void execute_kernels(struct PsiVector* state, const struct PsiKernel* kernels,
+						    size_t count, size_t num_qubits, struct PsiRuntimeConfig config)
+{
+	...
+}
+```
+
+### Includes
+
+Not sorted by the tool — order them yourself. Convention: system headers (`<...>`) first, a
+blank line, then project headers (`"..."`). `include/psi.h` orders its headers by dependency
+layer, not alphabetically.
+
+---
+
+## 3. Naming (enforced by clang-tidy)
+
+| Kind | Convention | Example |
+|------|-----------|---------|
+| Structs | `Psi` + PascalCase | `struct PsiQuantumRegister` |
+| Enums | `Psi` + PascalCase | `enum PsiGateType` |
+| Enum constants | `UPPER_CASE` | `PSI_GATE_H`, `PSI_SIMD_NEON` |
+| Functions | `snake_case` | `psi_add_complex`, `apply_kernel` |
+| Parameters / variables | `snake_case` | `num_qubits`, `target_bit` |
+| Macros / compile-time consts | `THIS_CASE` | `PSI_VERSION_MAJOR`, `INV_SQRT_2` |
+| Function-like macros | `psi_` + `snake_case` | `psi_matrix`, `psi_column_vector` |
+
+### Library prefix
+
+Every **public** symbol is namespaced:
+
+- functions → `psi_...`
+- structs → `Psi...`
+- macros / constants → `PSI_...`
+
+File-local `static` helpers are **not** prefixed (`apply_pair`, `build_pairs`, `op_to_kernel`).
+
+### Data-driven, action-first (not OOP)
+
+The operation leads and the data is an argument — think free functions over data, not methods
+bound to a type. `psi_<operation>_<type>(...)`:
+
+```c
+psi_add_complex(a, b)       /* not psi_complex_add   */
+psi_new_quantum_gate(...)   /* not gate_new(...)     */
+psi_apply_gate(&reg, ...)   /* verb first            */
+```
+
+*Not enforceable by tooling — maintained in review.*
+
+---
+
+## 4. Types
+
+- **`double` only** for floating point. No `float`, no fixed-width float aliases (their sizes
+  are not guaranteed across platforms).
+- **Fixed-width `<stdint.h>`** types (`uint32_t`, `int64_t`, …) over `int`/`unsigned`/`long`.
+- **`size_t`** for sizes, counts, and indices (the Rust `usize`).
+- **No `typedef` on structs, enums, or unions** — always spell `struct PsiComplex`,
+  `enum PsiGateType`.
+
+---
+
+## 5. Memory & ownership
+
+There is no garbage collector; ownership is explicit and follows a few rules.
+
+- Every type that owns a heap allocation has a matching `psi_new_*` / `psi_free_*` pair, and
+  `psi_free_*` takes a pointer and nulls the freed fields.
+
+```c
+struct PsiVector v = psi_new_vector(4, PSI_COLUMN_VECTOR);
+...
+psi_free_vector(&v);
+```
+
+- **The caller frees.** Functions that allocate and return a value (`psi_add_vector`,
+  `psi_clone_matrix`, the gate factories, the renderers' returned strings) transfer ownership
+  to the caller.
+- **Read-only args are passed by value**, which shares the underlying buffer and never frees
+  it. Mutators and destructors take a pointer.
+- **Constructors that take ownership** consume what they are handed (e.g. a gate takes
+  ownership of its matrix; `psi_apply_custom` takes ownership of the custom gate). Builders
+  that take an existing value they should not consume make a copy (`psi_new_vector_from`,
+  `psi_new_quantum_register_from`).
+- **Allocation idiom:** `sizeof *ptr`, not the repeated type.
+
+```c
+struct PsiKernel* out = malloc(count * sizeof *out);
+```
+
+- Verify with sanitizers during development:
+
+```bash
+clang -std=c17 -Iinclude -Isrc -Itester -fsanitize=address,undefined -g \
+	src/psi.c src/**/*.c tester/*.c -o /tmp/psi -lm && /tmp/psi
+```
+
+---
+
+## 6. Error handling
+
+- Preconditions and programmer errors (dimension mismatches, out-of-range indices) use
+  `assert`. The Rust original returned `Option`/panicked; in C these are `assert`s, since they
+  indicate caller bugs rather than recoverable conditions.
+- `malloc`/`calloc`/`realloc` results are `assert`ed non-NULL (with `|| size == 0` where a
+  zero-size allocation is legal).
+
+---
+
+## 7. Comments
+
+- **Avoid narration.** Code should read on its own; do not annotate ported steps or restate
+  what a line does.
+- **Short math-formula comments are welcome** where they clarify an expression:
+
+```c
+struct PsiComplex psi_mul_complex(struct PsiComplex a, struct PsiComplex b)
+{
+	// (a + bi)(c + di) = (ac - bd) + (ad + bc)i
+	return (struct PsiComplex){
+		a.real * b.real - a.imaginary * b.imaginary,
+		a.real * b.imaginary + a.imaginary * b.real,
+	};
+}
+```
+
+*Not enforceable by tooling — maintained in review.*
+
+---
+
+## 8. Project layout
+
+```
+include/          public API (mirrors src/ subtree)
+  psi.h           umbrella header — includes the whole public API
+  maths/          complex, vector, matrix, format, simd
+  core/           quantum_components, gates, custom_gate, classical_components,
+                  circuit, kernel, runtime, noise
+  visualizer/     renderer
+src/              implementation, same subtree; may also hold private headers
+tester/           assertion-based test suite
+```
+
+- Headers are `#pragma once` and included **root-relative**: `#include "maths/complex.h"`.
+- `include/psi.h` is the umbrella: `#include <psi.h>` pulls in the entire public API. Add each
+  new module's public header to it.
+- Private, implementation-only headers live under `src/` (e.g. `src/visualizer/grid.h`).
+
+---
+
+## 9. Testing
+
+- `tester/` is an assertion suite with subcommands
+  (`clifford`, `non-clifford`, `custom`, `kernels`, `simd`, `noise`, `all`, `help`).
+- Tests assert against known states (Bell, GHZ, rotations, fusion, noise purity, …) and check
+  that all runtimes agree. The process exit code reflects pass/fail (CI-friendly).
+
+```bash
+premake5 gmake && make config=release
+./bin/release-<system>/tester            # all suites
+./bin/release-<system>/tester noise      # one suite
+```
+
+---
+
+## Tooling summary
+
+| Concern | Tool | Command |
+|---------|------|---------|
+| Formatting | clang-format | `... | xargs clang-format --dry-run -Werror` |
+| Naming | clang-tidy | `clang-tidy <file> -- -std=c17 -Iinclude -Isrc` |
+| Build | premake5 + make | `premake5 gmake && make config=release` |
+| Tests | tester | `./bin/release-<system>/tester` |
+| Leaks / UB | sanitizers | build with `-fsanitize=address,undefined` |
diff --git a/include/core/circuit.h b/include/core/circuit.h
index bbc864c..7c0d371 100644
--- a/include/core/circuit.h
+++ b/include/core/circuit.h
@@ -46,14 +46,14 @@ struct PsiGateOp
 	double params[3];
 	size_t param_count;
 	size_t classical;
-	struct PsiCustomGate *custom;
-	size_t *custom_targets;
+	struct PsiCustomGate* custom;
+	size_t* custom_targets;
 	size_t custom_target_count;
 };
 
-const char *psi_gate_op_name(struct PsiGateOp op);
-const size_t *psi_gate_op_quantum_targets(const struct PsiGateOp *op, size_t *out_count);
-const size_t *psi_gate_op_classical_targets(const struct PsiGateOp *op, size_t *out_count);
+const char* psi_gate_op_name(struct PsiGateOp op);
+const size_t* psi_gate_op_quantum_targets(const struct PsiGateOp* op, size_t* out_count);
+const size_t* psi_gate_op_classical_targets(const struct PsiGateOp* op, size_t* out_count);
 bool psi_gate_op_is_measurement(struct PsiGateOp op);
 bool psi_gate_op_is_custom(struct PsiGateOp op);
 bool psi_gate_op_is_non_clifford(struct PsiGateOp op);
@@ -62,7 +62,7 @@ struct PsiQuantumCircuit
 {
 	size_t num_qubits;
 	size_t num_classical;
-	struct PsiGateOp *operations;
+	struct PsiGateOp* operations;
 	size_t operation_count;
 	size_t operation_capacity;
 	struct PsiVector computed_state;
@@ -70,41 +70,44 @@ struct PsiQuantumCircuit
 };
 
 struct PsiQuantumCircuit psi_new_quantum_circuit(size_t num_qubits);
-struct PsiQuantumCircuit psi_new_quantum_circuit_with_classical(size_t num_qubits, size_t num_classical);
-void psi_free_quantum_circuit(struct PsiQuantumCircuit *c);
-void psi_reset_circuit(struct PsiQuantumCircuit *c);
+struct PsiQuantumCircuit psi_new_quantum_circuit_with_classical(size_t num_qubits,
+                                                                size_t num_classical);
+void psi_free_quantum_circuit(struct PsiQuantumCircuit* c);
+void psi_reset_circuit(struct PsiQuantumCircuit* c);
 
-void psi_apply_h(struct PsiQuantumCircuit *c, size_t target);
-void psi_apply_x(struct PsiQuantumCircuit *c, size_t target);
-void psi_apply_y(struct PsiQuantumCircuit *c, size_t target);
-void psi_apply_z(struct PsiQuantumCircuit *c, size_t target);
-void psi_apply_s(struct PsiQuantumCircuit *c, size_t target);
-void psi_apply_t(struct PsiQuantumCircuit *c, size_t target);
-void psi_apply_sdg(struct PsiQuantumCircuit *c, size_t target);
-void psi_apply_tdg(struct PsiQuantumCircuit *c, size_t target);
-void psi_apply_sx(struct PsiQuantumCircuit *c, size_t target);
-void psi_apply_sxdg(struct PsiQuantumCircuit *c, size_t target);
+void psi_apply_h(struct PsiQuantumCircuit* c, size_t target);
+void psi_apply_x(struct PsiQuantumCircuit* c, size_t target);
+void psi_apply_y(struct PsiQuantumCircuit* c, size_t target);
+void psi_apply_z(struct PsiQuantumCircuit* c, size_t target);
+void psi_apply_s(struct PsiQuantumCircuit* c, size_t target);
+void psi_apply_t(struct PsiQuantumCircuit* c, size_t target);
+void psi_apply_sdg(struct PsiQuantumCircuit* c, size_t target);
+void psi_apply_tdg(struct PsiQuantumCircuit* c, size_t target);
+void psi_apply_sx(struct PsiQuantumCircuit* c, size_t target);
+void psi_apply_sxdg(struct PsiQuantumCircuit* c, size_t target);
 
-void psi_apply_rx(struct PsiQuantumCircuit *c, size_t target, double theta);
-void psi_apply_ry(struct PsiQuantumCircuit *c, size_t target, double theta);
-void psi_apply_rz(struct PsiQuantumCircuit *c, size_t target, double theta);
-void psi_apply_p(struct PsiQuantumCircuit *c, size_t target, double theta);
-void psi_apply_u1(struct PsiQuantumCircuit *c, size_t target, double lambda);
-void psi_apply_u2(struct PsiQuantumCircuit *c, size_t target, double phi, double lambda);
-void psi_apply_u3(struct PsiQuantumCircuit *c, size_t target, double theta, double phi, double lambda);
+void psi_apply_rx(struct PsiQuantumCircuit* c, size_t target, double theta);
+void psi_apply_ry(struct PsiQuantumCircuit* c, size_t target, double theta);
+void psi_apply_rz(struct PsiQuantumCircuit* c, size_t target, double theta);
+void psi_apply_p(struct PsiQuantumCircuit* c, size_t target, double theta);
+void psi_apply_u1(struct PsiQuantumCircuit* c, size_t target, double lambda);
+void psi_apply_u2(struct PsiQuantumCircuit* c, size_t target, double phi, double lambda);
+void psi_apply_u3(struct PsiQuantumCircuit* c, size_t target, double theta, double phi,
+                  double lambda);
 
-void psi_apply_cnot(struct PsiQuantumCircuit *c, size_t control, size_t target);
-void psi_apply_cz(struct PsiQuantumCircuit *c, size_t control, size_t target);
-void psi_apply_swap(struct PsiQuantumCircuit *c, size_t qubit1, size_t qubit2);
-void psi_apply_crx(struct PsiQuantumCircuit *c, size_t control, size_t target, double theta);
-void psi_apply_cry(struct PsiQuantumCircuit *c, size_t control, size_t target, double theta);
-void psi_apply_crz(struct PsiQuantumCircuit *c, size_t control, size_t target, double theta);
-void psi_apply_cp(struct PsiQuantumCircuit *c, size_t control, size_t target, double theta);
+void psi_apply_cnot(struct PsiQuantumCircuit* c, size_t control, size_t target);
+void psi_apply_cz(struct PsiQuantumCircuit* c, size_t control, size_t target);
+void psi_apply_swap(struct PsiQuantumCircuit* c, size_t qubit1, size_t qubit2);
+void psi_apply_crx(struct PsiQuantumCircuit* c, size_t control, size_t target, double theta);
+void psi_apply_cry(struct PsiQuantumCircuit* c, size_t control, size_t target, double theta);
+void psi_apply_crz(struct PsiQuantumCircuit* c, size_t control, size_t target, double theta);
+void psi_apply_cp(struct PsiQuantumCircuit* c, size_t control, size_t target, double theta);
 
-void psi_apply_ccnot(struct PsiQuantumCircuit *c, size_t control1, size_t control2, size_t target);
-void psi_apply_cswap(struct PsiQuantumCircuit *c, size_t control, size_t target1, size_t target2);
+void psi_apply_ccnot(struct PsiQuantumCircuit* c, size_t control1, size_t control2, size_t target);
+void psi_apply_cswap(struct PsiQuantumCircuit* c, size_t control, size_t target1, size_t target2);
 
-void psi_measure(struct PsiQuantumCircuit *c, size_t qubit, size_t classical);
-void psi_measure_all(struct PsiQuantumCircuit *c);
+void psi_measure(struct PsiQuantumCircuit* c, size_t qubit, size_t classical);
+void psi_measure_all(struct PsiQuantumCircuit* c);
 
-void psi_apply_custom(struct PsiQuantumCircuit *c, struct PsiCustomGate gate, const size_t *targets, size_t count);
+void psi_apply_custom(struct PsiQuantumCircuit* c, struct PsiCustomGate gate, const size_t* targets,
+                      size_t count);
diff --git a/include/core/classical_components.h b/include/core/classical_components.h
index 48714b4..53c00ef 100644
--- a/include/core/classical_components.h
+++ b/include/core/classical_components.h
@@ -5,17 +5,18 @@
 
 struct PsiClassicalBit
 {
-	const char *name;
+	const char* name;
 	bool state;
 };
 
 struct PsiClassicalRegister
 {
-	const char *name;
-	struct PsiClassicalBit *bits;
+	const char* name;
+	struct PsiClassicalBit* bits;
 	size_t num_bits;
 };
 
-struct PsiClassicalBit psi_new_classical_bit(const char *name, bool state);
-struct PsiClassicalRegister psi_new_classical_register(const char *name, const char **names, size_t count);
-void psi_free_classical_register(struct PsiClassicalRegister *reg);
+struct PsiClassicalBit psi_new_classical_bit(const char* name, bool state);
+struct PsiClassicalRegister psi_new_classical_register(const char* name, const char** names,
+                                                       size_t count);
+void psi_free_classical_register(struct PsiClassicalRegister* reg);
diff --git a/include/core/custom_gate.h b/include/core/custom_gate.h
index bc22956..efda36f 100644
--- a/include/core/custom_gate.h
+++ b/include/core/custom_gate.h
@@ -35,7 +35,7 @@ enum PsiCustomGateKind
 
 struct PsiCustomGate
 {
-	const char *name;
+	const char* name;
 	size_t num_qubits;
 	enum PsiCustomGateKind kind;
 	union
@@ -43,14 +43,16 @@ struct PsiCustomGate
 		struct PsiMatrix matrix;
 		struct
 		{
-			struct PsiCompositeGateOp *ops;
+			struct PsiCompositeGateOp* ops;
 			size_t count;
 		} composite;
 	} definition;
 };
 
-struct PsiCustomGate psi_new_custom_gate_from_matrix(const char *name, struct PsiMatrix matrix);
-struct PsiCustomGate psi_new_custom_gate_from_composite(const char *name, size_t num_qubits, const struct PsiCompositeGateOp *ops, size_t op_count);
-void psi_free_custom_gate(struct PsiCustomGate *gate);
+struct PsiCustomGate psi_new_custom_gate_from_matrix(const char* name, struct PsiMatrix matrix);
+struct PsiCustomGate psi_new_custom_gate_from_composite(const char* name, size_t num_qubits,
+                                                        const struct PsiCompositeGateOp* ops,
+                                                        size_t op_count);
+void psi_free_custom_gate(struct PsiCustomGate* gate);
 
 struct PsiQuantumGate psi_to_quantum_gate(struct PsiCustomGate gate);
diff --git a/include/core/kernel.h b/include/core/kernel.h
index d7976c3..b8b19ea 100644
--- a/include/core/kernel.h
+++ b/include/core/kernel.h
@@ -16,39 +16,40 @@ enum PsiGateType
 struct PsiKernel
 {
 	struct PsiMatrix matrix;
-	size_t *targets;
+	size_t* targets;
 	size_t target_count;
-	char *name;
+	char* name;
 	enum PsiGateType gate_type;
 };
 
-struct PsiKernel psi_new_kernel(const char *name, struct PsiMatrix matrix, const size_t *targets, size_t target_count);
+struct PsiKernel psi_new_kernel(const char* name, struct PsiMatrix matrix, const size_t* targets,
+                                size_t target_count);
 struct PsiKernel psi_clone_kernel(struct PsiKernel kernel);
-void psi_free_kernel(struct PsiKernel *kernel);
+void psi_free_kernel(struct PsiKernel* kernel);
 
 bool psi_kernels_share_qubits(struct PsiKernel a, struct PsiKernel b);
 bool psi_kernels_commute(struct PsiKernel a, struct PsiKernel b);
 bool psi_kernels_can_fuse(struct PsiKernel a, struct PsiKernel b);
-bool psi_fuse_kernels(struct PsiKernel a, struct PsiKernel b, struct PsiKernel *out);
+bool psi_fuse_kernels(struct PsiKernel a, struct PsiKernel b, struct PsiKernel* out);
 
 struct PsiKernelBatch
 {
-	struct PsiKernel *kernels;
+	struct PsiKernel* kernels;
 	size_t count;
 	size_t capacity;
 	size_t num_qubits;
 };
 
 struct PsiKernelBatch psi_new_kernel_batch(size_t num_qubits);
-void psi_free_kernel_batch(struct PsiKernelBatch *batch);
-void psi_add_kernel(struct PsiKernelBatch *batch, struct PsiKernel kernel);
-void psi_optimize_kernel_batch(struct PsiKernelBatch *batch);
-void psi_execute_kernel_batch(struct PsiKernelBatch batch, struct PsiVector *state);
-void psi_apply_kernel(struct PsiVector *state, struct PsiKernel kernel, size_t num_qubits);
+void psi_free_kernel_batch(struct PsiKernelBatch* batch);
+void psi_add_kernel(struct PsiKernelBatch* batch, struct PsiKernel kernel);
+void psi_optimize_kernel_batch(struct PsiKernelBatch* batch);
+void psi_execute_kernel_batch(struct PsiKernelBatch batch, struct PsiVector* state);
+void psi_apply_kernel(struct PsiVector* state, struct PsiKernel kernel, size_t num_qubits);
 
 struct PsiExecutionLayer
 {
-	struct PsiKernel *kernels;
+	struct PsiKernel* kernels;
 	size_t count;
 	size_t capacity;
 };
@@ -65,10 +66,10 @@ struct PsiKernelStats
 
 struct PsiStructureAwareBatch
 {
-	struct PsiKernel *kernels;
+	struct PsiKernel* kernels;
 	size_t count;
 	size_t capacity;
-	struct PsiExecutionLayer *layers;
+	struct PsiExecutionLayer* layers;
 	size_t layer_count;
 	size_t layer_capacity;
 	size_t num_qubits;
@@ -76,9 +77,11 @@ struct PsiStructureAwareBatch
 };
 
 struct PsiStructureAwareBatch psi_new_structure_aware_batch(size_t num_qubits);
-void psi_free_structure_aware_batch(struct PsiStructureAwareBatch *batch);
-void psi_add_structure_aware_kernel(struct PsiStructureAwareBatch *batch, struct PsiKernel kernel);
-void psi_optimize_structure_aware_batch(struct PsiStructureAwareBatch *batch);
-void psi_execute_structure_aware_batch(struct PsiStructureAwareBatch batch, struct PsiVector *state);
-void psi_execute_structure_aware_batch_layered(struct PsiStructureAwareBatch batch, struct PsiVector *state);
+void psi_free_structure_aware_batch(struct PsiStructureAwareBatch* batch);
+void psi_add_structure_aware_kernel(struct PsiStructureAwareBatch* batch, struct PsiKernel kernel);
+void psi_optimize_structure_aware_batch(struct PsiStructureAwareBatch* batch);
+void psi_execute_structure_aware_batch(struct PsiStructureAwareBatch batch,
+                                       struct PsiVector* state);
+void psi_execute_structure_aware_batch_layered(struct PsiStructureAwareBatch batch,
+                                               struct PsiVector* state);
 struct PsiKernelStats psi_structure_aware_batch_stats(struct PsiStructureAwareBatch batch);
diff --git a/include/core/noise.h b/include/core/noise.h
index 6232d8e..320a4ab 100644
--- a/include/core/noise.h
+++ b/include/core/noise.h
@@ -9,22 +9,24 @@
 struct PsiKrausOperator
 {
 	struct PsiMatrix matrix;
-	const char *name;
+	const char* name;
 };
 
-struct PsiKrausOperator psi_new_kraus_operator(const char *name, struct PsiMatrix matrix);
-void psi_free_kraus_operator(struct PsiKrausOperator *op);
+struct PsiKrausOperator psi_new_kraus_operator(const char* name, struct PsiMatrix matrix);
+void psi_free_kraus_operator(struct PsiKrausOperator* op);
 
 struct PsiNoiseChannel
 {
-	const char *name;
-	struct PsiKrausOperator *operators;
+	const char* name;
+	struct PsiKrausOperator* operators;
 	size_t operator_count;
 	size_t num_qubits;
 };
 
-struct PsiNoiseChannel psi_new_noise_channel(const char *name, const struct PsiKrausOperator *operators, size_t count, size_t num_qubits);
-void psi_free_noise_channel(struct PsiNoiseChannel *channel);
+struct PsiNoiseChannel psi_new_noise_channel(const char* name,
+                                             const struct PsiKrausOperator* operators, size_t count,
+                                             size_t num_qubits);
+void psi_free_noise_channel(struct PsiNoiseChannel* channel);
 
 struct PsiNoiseChannel psi_depolarising_channel(double p);
 struct PsiNoiseChannel psi_amplitude_damping_channel(double gamma);
@@ -36,22 +38,27 @@ struct PsiNoiseChannel psi_generalised_amplitude_damping_channel(double p, doubl
 
 struct PsiDensityMatrix
 {
-	struct PsiComplex *data;
+	struct PsiComplex* data;
 	size_t dim;
 	size_t num_qubits;
 };
 
 struct PsiDensityMatrix psi_new_density_matrix(size_t num_qubits);
-struct PsiDensityMatrix psi_new_density_matrix_from_state(const struct PsiComplex *state, size_t len);
-void psi_free_density_matrix(struct PsiDensityMatrix *dm);
+struct PsiDensityMatrix psi_new_density_matrix_from_state(const struct PsiComplex* state,
+                                                          size_t len);
+void psi_free_density_matrix(struct PsiDensityMatrix* dm);
 
 struct PsiComplex psi_get_density_matrix(struct PsiDensityMatrix dm, size_t row, size_t col);
-void psi_set_density_matrix(struct PsiDensityMatrix *dm, size_t row, size_t col, struct PsiComplex value);
+void psi_set_density_matrix(struct PsiDensityMatrix* dm, size_t row, size_t col,
+                            struct PsiComplex value);
 struct PsiComplex psi_trace_density_matrix(struct PsiDensityMatrix dm);
 double psi_purity_density_matrix(struct PsiDensityMatrix dm);
 bool psi_is_pure_density_matrix(struct PsiDensityMatrix dm, double tolerance);
-void psi_density_matrix_probabilities(struct PsiDensityMatrix dm, double *out);
-void psi_apply_unitary_density_matrix(struct PsiDensityMatrix *dm, struct PsiMatrix gate, const size_t *targets, size_t target_count);
-void psi_apply_noise_channel(struct PsiDensityMatrix *dm, struct PsiNoiseChannel channel, size_t target);
-double psi_measure_probability_density_matrix(struct PsiDensityMatrix dm, size_t qubit, size_t outcome);
-double psi_fidelity_density_matrix(struct PsiDensityMatrix dm, const struct PsiComplex *state);
+void psi_density_matrix_probabilities(struct PsiDensityMatrix dm, double* out);
+void psi_apply_unitary_density_matrix(struct PsiDensityMatrix* dm, struct PsiMatrix gate,
+                                      const size_t* targets, size_t target_count);
+void psi_apply_noise_channel(struct PsiDensityMatrix* dm, struct PsiNoiseChannel channel,
+                             size_t target);
+double psi_measure_probability_density_matrix(struct PsiDensityMatrix dm, size_t qubit,
+                                              size_t outcome);
+double psi_fidelity_density_matrix(struct PsiDensityMatrix dm, const struct PsiComplex* state);
diff --git a/include/core/quantum_components.h b/include/core/quantum_components.h
index b238b5a..63797cf 100644
--- a/include/core/quantum_components.h
+++ b/include/core/quantum_components.h
@@ -10,34 +10,38 @@ struct PsiVector psi_new_state_1(void);
 
 struct PsiQuantumGate
 {
-	const char *name;
+	const char* name;
 	struct PsiMatrix matrix;
 	size_t num_qubits;
 };
 
-struct PsiQuantumGate psi_new_quantum_gate(const char *name, struct PsiMatrix matrix, size_t num_qubits);
-struct PsiQuantumGate psi_new_quantum_gate_from_matrix(const char *name, struct PsiMatrix matrix);
-void psi_free_quantum_gate(struct PsiQuantumGate *gate);
+struct PsiQuantumGate psi_new_quantum_gate(const char* name, struct PsiMatrix matrix,
+                                           size_t num_qubits);
+struct PsiQuantumGate psi_new_quantum_gate_from_matrix(const char* name, struct PsiMatrix matrix);
+void psi_free_quantum_gate(struct PsiQuantumGate* gate);
 
 struct PsiQuantumBit
 {
-	const char *name;
+	const char* name;
 	struct PsiVector state;
 };
 
-struct PsiQuantumBit psi_new_quantum_bit(const char *name, struct PsiVector state);
-void psi_free_quantum_bit(struct PsiQuantumBit *bit);
+struct PsiQuantumBit psi_new_quantum_bit(const char* name, struct PsiVector state);
+void psi_free_quantum_bit(struct PsiQuantumBit* bit);
 
 struct PsiQuantumRegister
 {
-	const char *name;
+	const char* name;
 	struct PsiVector state_vector;
-	struct PsiQuantumBit *qubits;
+	struct PsiQuantumBit* qubits;
 	size_t num_qubits;
 };
 
-struct PsiQuantumRegister psi_new_quantum_register(const char *name, const char **names, size_t count);
-struct PsiQuantumRegister psi_new_quantum_register_from(const char *name, const struct PsiQuantumBit *bits, size_t count);
-void psi_free_quantum_register(struct PsiQuantumRegister *reg);
+struct PsiQuantumRegister psi_new_quantum_register(const char* name, const char** names,
+                                                   size_t count);
+struct PsiQuantumRegister
+psi_new_quantum_register_from(const char* name, const struct PsiQuantumBit* bits, size_t count);
+void psi_free_quantum_register(struct PsiQuantumRegister* reg);
 
-void psi_apply_gate(struct PsiQuantumRegister *reg, struct PsiQuantumGate gate, const size_t *targets, size_t target_count);
+void psi_apply_gate(struct PsiQuantumRegister* reg, struct PsiQuantumGate gate,
+                    const size_t* targets, size_t target_count);
diff --git a/include/core/runtime.h b/include/core/runtime.h
index 9e7427e..2a6f739 100644
--- a/include/core/runtime.h
+++ b/include/core/runtime.h
@@ -33,11 +33,15 @@ struct PsiRuntimeConfig psi_new_runtime_config(void);
 struct PsiRuntimeConfig psi_optimal_runtime_config(void);
 struct PsiRuntimeConfig psi_runtime_to_config(enum PsiRuntime runtime);
 
-struct PsiVector psi_compute_runtime_config(struct PsiRuntimeConfig config, size_t num_qubits, const struct PsiGateOp *operations, size_t op_count);
-struct PsiVector psi_compute_runtime(enum PsiRuntime runtime, size_t num_qubits, const struct PsiGateOp *operations, size_t op_count);
+struct PsiVector psi_compute_runtime_config(struct PsiRuntimeConfig config, size_t num_qubits,
+                                            const struct PsiGateOp* operations, size_t op_count);
+struct PsiVector psi_compute_runtime(enum PsiRuntime runtime, size_t num_qubits,
+                                     const struct PsiGateOp* operations, size_t op_count);
 
-const struct PsiVector *psi_compute_circuit(struct PsiQuantumCircuit *circuit);
-const struct PsiVector *psi_compute_circuit_with(struct PsiQuantumCircuit *circuit, enum PsiRuntime runtime);
-const struct PsiVector *psi_compute_circuit_with_config(struct PsiQuantumCircuit *circuit, struct PsiRuntimeConfig config);
-double psi_circuit_probability(struct PsiQuantumCircuit *circuit, size_t state_index);
-void psi_circuit_probabilities(struct PsiQuantumCircuit *circuit, double *out);
+const struct PsiVector* psi_compute_circuit(struct PsiQuantumCircuit* circuit);
+const struct PsiVector* psi_compute_circuit_with(struct PsiQuantumCircuit* circuit,
+                                                 enum PsiRuntime runtime);
+const struct PsiVector* psi_compute_circuit_with_config(struct PsiQuantumCircuit* circuit,
+                                                        struct PsiRuntimeConfig config);
+double psi_circuit_probability(struct PsiQuantumCircuit* circuit, size_t state_index);
+void psi_circuit_probabilities(struct PsiQuantumCircuit* circuit, double* out);
diff --git a/include/maths/format.h b/include/maths/format.h
index 65c1680..a11e3be 100644
--- a/include/maths/format.h
+++ b/include/maths/format.h
@@ -4,5 +4,5 @@
 
 #include "maths/complex.h"
 
-char *psi_format_amplitude(struct PsiComplex c, char *out, size_t cap);
-char *psi_format_probability(double p, char *out, size_t cap);
+char* psi_format_amplitude(struct PsiComplex c, char* out, size_t cap);
+char* psi_format_probability(double p, char* out, size_t cap);
diff --git a/include/maths/matrix.h b/include/maths/matrix.h
index c4e5530..ff996d6 100644
--- a/include/maths/matrix.h
+++ b/include/maths/matrix.h
@@ -7,22 +7,22 @@
 
 struct PsiMatrix
 {
-	struct PsiComplex *data;
+	struct PsiComplex* data;
 	size_t rows;
 	size_t cols;
 };
 
 struct PsiMatrix psi_new_matrix(size_t rows, size_t cols);
-struct PsiMatrix psi_new_matrix_from(const struct PsiComplex *data, size_t rows, size_t cols);
+struct PsiMatrix psi_new_matrix_from(const struct PsiComplex* data, size_t rows, size_t cols);
 struct PsiMatrix psi_identity_matrix(size_t size);
 struct PsiMatrix psi_clone_matrix(struct PsiMatrix m);
-void psi_free_matrix(struct PsiMatrix *m);
+void psi_free_matrix(struct PsiMatrix* m);
 
-#define psi_matrix(rows, cols, ...) \
+#define psi_matrix(rows, cols, ...)                                                                \
 	psi_new_matrix_from((struct PsiComplex[]){ __VA_ARGS__ }, rows, cols)
 
 struct PsiComplex psi_get_matrix(struct PsiMatrix m, size_t row, size_t col);
-void psi_set_matrix(struct PsiMatrix *m, size_t row, size_t col, struct PsiComplex value);
+void psi_set_matrix(struct PsiMatrix* m, size_t row, size_t col, struct PsiComplex value);
 
 struct PsiMatrix psi_dot_matrix(struct PsiMatrix a, struct PsiMatrix b);
 struct PsiMatrix psi_kronecker_matrix(struct PsiMatrix a, struct PsiMatrix b);
diff --git a/include/maths/simd.h b/include/maths/simd.h
new file mode 100644
index 0000000..38d0171
--- /dev/null
+++ b/include/maths/simd.h
@@ -0,0 +1,22 @@
+#pragma once
+
+#include <stddef.h>
+
+#include "maths/complex.h"
+
+enum PsiSimdCapability
+{
+	PSI_SIMD_NONE,
+	PSI_SIMD_AVX2,
+	PSI_SIMD_AVX512,
+	PSI_SIMD_NEON,
+};
+
+enum PsiSimdCapability psi_detect_simd(void);
+const char* psi_simd_name(enum PsiSimdCapability cap);
+
+void psi_apply_single_qubit_gate_simd(struct PsiComplex* state, const struct PsiComplex gate[2][2],
+                                      size_t target, size_t num_qubits);
+void psi_apply_single_qubit_gate_simd_parallel(struct PsiComplex* state,
+                                               const struct PsiComplex gate[2][2], size_t target,
+                                               size_t num_qubits);
diff --git a/include/maths/vector.h b/include/maths/vector.h
index d452ef1..0136715 100644
--- a/include/maths/vector.h
+++ b/include/maths/vector.h
@@ -12,27 +12,28 @@ enum PsiVectorKind
 
 struct PsiVector
 {
-	struct PsiComplex *data;
+	struct PsiComplex* data;
 	size_t size;
 	enum PsiVectorKind kind;
 };
 
 struct PsiVector psi_new_vector(size_t size, enum PsiVectorKind kind);
-struct PsiVector psi_new_vector_from(const struct PsiComplex *data, size_t size, enum PsiVectorKind kind);
+struct PsiVector psi_new_vector_from(const struct PsiComplex* data, size_t size,
+                                     enum PsiVectorKind kind);
 struct PsiVector psi_clone_vector(struct PsiVector v);
-void psi_free_vector(struct PsiVector *v);
+void psi_free_vector(struct PsiVector* v);
 
-#define psi_row_vector(...) \
-	psi_new_vector_from((struct PsiComplex[]){ __VA_ARGS__ }, \
-		sizeof((struct PsiComplex[]){ __VA_ARGS__ }) / sizeof(struct PsiComplex), \
-		PSI_ROW_VECTOR)
-#define psi_column_vector(...) \
-	psi_new_vector_from((struct PsiComplex[]){ __VA_ARGS__ }, \
-		sizeof((struct PsiComplex[]){ __VA_ARGS__ }) / sizeof(struct PsiComplex), \
-		PSI_COLUMN_VECTOR)
+#define psi_row_vector(...)                                                                        \
+	psi_new_vector_from((struct PsiComplex[]){ __VA_ARGS__ },                                      \
+	                    sizeof((struct PsiComplex[]){ __VA_ARGS__ }) / sizeof(struct PsiComplex),  \
+	                    PSI_ROW_VECTOR)
+#define psi_column_vector(...)                                                                     \
+	psi_new_vector_from((struct PsiComplex[]){ __VA_ARGS__ },                                      \
+	                    sizeof((struct PsiComplex[]){ __VA_ARGS__ }) / sizeof(struct PsiComplex),  \
+	                    PSI_COLUMN_VECTOR)
 
 struct PsiComplex psi_get_vector(struct PsiVector v, size_t index);
-void psi_set_vector(struct PsiVector *v, size_t index, struct PsiComplex value);
+void psi_set_vector(struct PsiVector* v, size_t index, struct PsiComplex value);
 
 struct PsiComplex psi_dot_vector(struct PsiVector a, struct PsiVector b);
 struct PsiComplex psi_norm_vector(struct PsiVector v);
diff --git a/include/psi.h b/include/psi.h
index 6bc673a..afd7d89 100644
--- a/include/psi.h
+++ b/include/psi.h
@@ -4,12 +4,13 @@
 #define PSI_VERSION_MINOR 1
 #define PSI_VERSION_PATCH 0
 
-const char *psi_version(void);
+const char* psi_version(void);
 
 #include "maths/complex.h"
 #include "maths/vector.h"
 #include "maths/matrix.h"
 #include "maths/format.h"
+#include "maths/simd.h"
 
 #include "core/quantum_components.h"
 #include "core/gates.h"
diff --git a/include/visualizer/renderer.h b/include/visualizer/renderer.h
index faabf2c..c3b1a1c 100644
--- a/include/visualizer/renderer.h
+++ b/include/visualizer/renderer.h
@@ -2,5 +2,5 @@
 
 #include "core/circuit.h"
 
-char *psi_render_circuit_vertical(const struct PsiQuantumCircuit *circuit);
-char *psi_render_circuit_horizontal(const struct PsiQuantumCircuit *circuit);
+char* psi_render_circuit_vertical(const struct PsiQuantumCircuit* circuit);
+char* psi_render_circuit_horizontal(const struct PsiQuantumCircuit* circuit);
diff --git a/src/core/circuit.c b/src/core/circuit.c
index 1869ee6..8009db5 100644
--- a/src/core/circuit.c
+++ b/src/core/circuit.c
@@ -4,7 +4,7 @@
 #include <stdlib.h>
 #include <string.h>
 
-const char *psi_gate_op_name(struct PsiGateOp op)
+const char* psi_gate_op_name(struct PsiGateOp op)
 {
 	switch (op.kind)
 	{
@@ -41,7 +41,7 @@ const char *psi_gate_op_name(struct PsiGateOp op)
 	return "?";
 }
 
-const size_t *psi_gate_op_quantum_targets(const struct PsiGateOp *op, size_t *out_count)
+const size_t* psi_gate_op_quantum_targets(const struct PsiGateOp* op, size_t* out_count)
 {
 	if (op->kind == PSI_GATE_CUSTOM)
 	{
@@ -53,7 +53,7 @@ const size_t *psi_gate_op_quantum_targets(const struct PsiGateOp *op, size_t *ou
 	return op->qubits;
 }
 
-const size_t *psi_gate_op_classical_targets(const struct PsiGateOp *op, size_t *out_count)
+const size_t* psi_gate_op_classical_targets(const struct PsiGateOp* op, size_t* out_count)
 {
 	if (op->kind == PSI_GATE_MEASURE)
 	{
@@ -93,10 +93,8 @@ bool psi_gate_op_is_non_clifford(struct PsiGateOp op)
 		case PSI_GATE_CRX:
 		case PSI_GATE_CRY:
 		case PSI_GATE_CRZ:
-		case PSI_GATE_CP:
-			return true;
-		default:
-			return false;
+		case PSI_GATE_CP: return true;
+		default: return false;
 	}
 }
 
@@ -105,7 +103,8 @@ struct PsiQuantumCircuit psi_new_quantum_circuit(size_t num_qubits)
 	return psi_new_quantum_circuit_with_classical(num_qubits, 0);
 }
 
-struct PsiQuantumCircuit psi_new_quantum_circuit_with_classical(size_t num_qubits, size_t num_classical)
+struct PsiQuantumCircuit psi_new_quantum_circuit_with_classical(size_t num_qubits,
+                                                                size_t num_classical)
 {
 	struct PsiQuantumCircuit c;
 	c.num_qubits = num_qubits;
@@ -119,11 +118,11 @@ struct PsiQuantumCircuit psi_new_quantum_circuit_with_classical(size_t num_qubit
 	return c;
 }
 
-static void free_operations(struct PsiQuantumCircuit *c)
+static void free_operations(struct PsiQuantumCircuit* c)
 {
 	for (size_t i = 0; i < c->operation_count; i++)
 	{
-		struct PsiGateOp *op = &c->operations[i];
+		struct PsiGateOp* op = &c->operations[i];
 		if (op->kind != PSI_GATE_CUSTOM)
 			continue;
 
@@ -133,7 +132,7 @@ static void free_operations(struct PsiQuantumCircuit *c)
 	}
 }
 
-void psi_free_quantum_circuit(struct PsiQuantumCircuit *c)
+void psi_free_quantum_circuit(struct PsiQuantumCircuit* c)
 {
 	free_operations(c);
 	free(c->operations);
@@ -144,7 +143,7 @@ void psi_free_quantum_circuit(struct PsiQuantumCircuit *c)
 	c->is_computed = false;
 }
 
-void psi_reset_circuit(struct PsiQuantumCircuit *c)
+void psi_reset_circuit(struct PsiQuantumCircuit* c)
 {
 	free_operations(c);
 	c->operation_count = 0;
@@ -153,7 +152,7 @@ void psi_reset_circuit(struct PsiQuantumCircuit *c)
 	c->is_computed = false;
 }
 
-static struct PsiGateOp *append_op(struct PsiQuantumCircuit *c)
+static struct PsiGateOp* append_op(struct PsiQuantumCircuit* c)
 {
 	if (c->operation_count == c->operation_capacity)
 	{
@@ -163,24 +162,25 @@ static struct PsiGateOp *append_op(struct PsiQuantumCircuit *c)
 		c->operation_capacity = new_capacity;
 	}
 
-	struct PsiGateOp *op = &c->operations[c->operation_count++];
+	struct PsiGateOp* op = &c->operations[c->operation_count++];
 	memset(op, 0, sizeof(*op));
 	c->is_computed = false;
 
 	return op;
 }
 
-static void push_1q(struct PsiQuantumCircuit *c, enum PsiGateKind kind, size_t target)
+static void push_1q(struct PsiQuantumCircuit* c, enum PsiGateKind kind, size_t target)
 {
-	struct PsiGateOp *op = append_op(c);
+	struct PsiGateOp* op = append_op(c);
 	op->kind = kind;
 	op->qubits[0] = target;
 	op->qubit_count = 1;
 }
 
-static void push_1q_1p(struct PsiQuantumCircuit *c, enum PsiGateKind kind, size_t target, double theta)
+static void push_1q_1p(struct PsiQuantumCircuit* c, enum PsiGateKind kind, size_t target,
+                       double theta)
 {
-	struct PsiGateOp *op = append_op(c);
+	struct PsiGateOp* op = append_op(c);
 	op->kind = kind;
 	op->qubits[0] = target;
 	op->qubit_count = 1;
@@ -188,18 +188,19 @@ static void push_1q_1p(struct PsiQuantumCircuit *c, enum PsiGateKind kind, size_
 	op->param_count = 1;
 }
 
-static void push_2q(struct PsiQuantumCircuit *c, enum PsiGateKind kind, size_t a, size_t b)
+static void push_2q(struct PsiQuantumCircuit* c, enum PsiGateKind kind, size_t a, size_t b)
 {
-	struct PsiGateOp *op = append_op(c);
+	struct PsiGateOp* op = append_op(c);
 	op->kind = kind;
 	op->qubits[0] = a;
 	op->qubits[1] = b;
 	op->qubit_count = 2;
 }
 
-static void push_2q_1p(struct PsiQuantumCircuit *c, enum PsiGateKind kind, size_t control, size_t target, double theta)
+static void push_2q_1p(struct PsiQuantumCircuit* c, enum PsiGateKind kind, size_t control,
+                       size_t target, double theta)
 {
-	struct PsiGateOp *op = append_op(c);
+	struct PsiGateOp* op = append_op(c);
 	op->kind = kind;
 	op->qubits[0] = control;
 	op->qubits[1] = target;
@@ -208,9 +209,10 @@ static void push_2q_1p(struct PsiQuantumCircuit *c, enum PsiGateKind kind, size_
 	op->param_count = 1;
 }
 
-static void push_3q(struct PsiQuantumCircuit *c, enum PsiGateKind kind, size_t a, size_t b, size_t d)
+static void push_3q(struct PsiQuantumCircuit* c, enum PsiGateKind kind, size_t a, size_t b,
+                    size_t d)
 {
-	struct PsiGateOp *op = append_op(c);
+	struct PsiGateOp* op = append_op(c);
 	op->kind = kind;
 	op->qubits[0] = a;
 	op->qubits[1] = b;
@@ -218,26 +220,71 @@ static void push_3q(struct PsiQuantumCircuit *c, enum PsiGateKind kind, size_t a
 	op->qubit_count = 3;
 }
 
-void psi_apply_h(struct PsiQuantumCircuit *c, size_t target) { push_1q(c, PSI_GATE_H, target); }
-void psi_apply_x(struct PsiQuantumCircuit *c, size_t target) { push_1q(c, PSI_GATE_X, target); }
-void psi_apply_y(struct PsiQuantumCircuit *c, size_t target) { push_1q(c, PSI_GATE_Y, target); }
-void psi_apply_z(struct PsiQuantumCircuit *c, size_t target) { push_1q(c, PSI_GATE_Z, target); }
-void psi_apply_s(struct PsiQuantumCircuit *c, size_t target) { push_1q(c, PSI_GATE_S, target); }
-void psi_apply_t(struct PsiQuantumCircuit *c, size_t target) { push_1q(c, PSI_GATE_T, target); }
-void psi_apply_sdg(struct PsiQuantumCircuit *c, size_t target) { push_1q(c, PSI_GATE_SDG, target); }
-void psi_apply_tdg(struct PsiQuantumCircuit *c, size_t target) { push_1q(c, PSI_GATE_TDG, target); }
-void psi_apply_sx(struct PsiQuantumCircuit *c, size_t target) { push_1q(c, PSI_GATE_SX, target); }
-void psi_apply_sxdg(struct PsiQuantumCircuit *c, size_t target) { push_1q(c, PSI_GATE_SXDG, target); }
+void psi_apply_h(struct PsiQuantumCircuit* c, size_t target)
+{
+	push_1q(c, PSI_GATE_H, target);
+}
+void psi_apply_x(struct PsiQuantumCircuit* c, size_t target)
+{
+	push_1q(c, PSI_GATE_X, target);
+}
+void psi_apply_y(struct PsiQuantumCircuit* c, size_t target)
+{
+	push_1q(c, PSI_GATE_Y, target);
+}
+void psi_apply_z(struct PsiQuantumCircuit* c, size_t target)
+{
+	push_1q(c, PSI_GATE_Z, target);
+}
+void psi_apply_s(struct PsiQuantumCircuit* c, size_t target)
+{
+	push_1q(c, PSI_GATE_S, target);
+}
+void psi_apply_t(struct PsiQuantumCircuit* c, size_t target)
+{
+	push_1q(c, PSI_GATE_T, target);
+}
+void psi_apply_sdg(struct PsiQuantumCircuit* c, size_t target)
+{
+	push_1q(c, PSI_GATE_SDG, target);
+}
+void psi_apply_tdg(struct PsiQuantumCircuit* c, size_t target)
+{
+	push_1q(c, PSI_GATE_TDG, target);
+}
+void psi_apply_sx(struct PsiQuantumCircuit* c, size_t target)
+{
+	push_1q(c, PSI_GATE_SX, target);
+}
+void psi_apply_sxdg(struct PsiQuantumCircuit* c, size_t target)
+{
+	push_1q(c, PSI_GATE_SXDG, target);
+}
 
-void psi_apply_rx(struct PsiQuantumCircuit *c, size_t target, double theta) { push_1q_1p(c, PSI_GATE_RX, target, theta); }
-void psi_apply_ry(struct PsiQuantumCircuit *c, size_t target, double theta) { push_1q_1p(c, PSI_GATE_RY, target, theta); }
-void psi_apply_rz(struct PsiQuantumCircuit *c, size_t target, double theta) { push_1q_1p(c, PSI_GATE_RZ, target, theta); }
-void psi_apply_p(struct PsiQuantumCircuit *c, size_t target, double theta) { push_1q_1p(c, PSI_GATE_P, target, theta); }
-void psi_apply_u1(struct PsiQuantumCircuit *c, size_t target, double lambda) { push_1q_1p(c, PSI_GATE_U1, target, lambda); }
+void psi_apply_rx(struct PsiQuantumCircuit* c, size_t target, double theta)
+{
+	push_1q_1p(c, PSI_GATE_RX, target, theta);
+}
+void psi_apply_ry(struct PsiQuantumCircuit* c, size_t target, double theta)
+{
+	push_1q_1p(c, PSI_GATE_RY, target, theta);
+}
+void psi_apply_rz(struct PsiQuantumCircuit* c, size_t target, double theta)
+{
+	push_1q_1p(c, PSI_GATE_RZ, target, theta);
+}
+void psi_apply_p(struct PsiQuantumCircuit* c, size_t target, double theta)
+{
+	push_1q_1p(c, PSI_GATE_P, target, theta);
+}
+void psi_apply_u1(struct PsiQuantumCircuit* c, size_t target, double lambda)
+{
+	push_1q_1p(c, PSI_GATE_U1, target, lambda);
+}
 
-void psi_apply_u2(struct PsiQuantumCircuit *c, size_t target, double phi, double lambda)
+void psi_apply_u2(struct PsiQuantumCircuit* c, size_t target, double phi, double lambda)
 {
-	struct PsiGateOp *op = append_op(c);
+	struct PsiGateOp* op = append_op(c);
 	op->kind = PSI_GATE_U2;
 	op->qubits[0] = target;
 	op->qubit_count = 1;
@@ -246,9 +293,10 @@ void psi_apply_u2(struct PsiQuantumCircuit *c, size_t target, double phi, double
 	op->param_count = 2;
 }
 
-void psi_apply_u3(struct PsiQuantumCircuit *c, size_t target, double theta, double phi, double lambda)
+void psi_apply_u3(struct PsiQuantumCircuit* c, size_t target, double theta, double phi,
+                  double lambda)
 {
-	struct PsiGateOp *op = append_op(c);
+	struct PsiGateOp* op = append_op(c);
 	op->kind = PSI_GATE_U3;
 	op->qubits[0] = target;
 	op->qubit_count = 1;
@@ -258,47 +306,75 @@ void psi_apply_u3(struct PsiQuantumCircuit *c, size_t target, double theta, doub
 	op->param_count = 3;
 }
 
-void psi_apply_cnot(struct PsiQuantumCircuit *c, size_t control, size_t target) { push_2q(c, PSI_GATE_CNOT, control, target); }
-void psi_apply_cz(struct PsiQuantumCircuit *c, size_t control, size_t target) { push_2q(c, PSI_GATE_CZ, control, target); }
-void psi_apply_swap(struct PsiQuantumCircuit *c, size_t qubit1, size_t qubit2) { push_2q(c, PSI_GATE_SWAP, qubit1, qubit2); }
-void psi_apply_crx(struct PsiQuantumCircuit *c, size_t control, size_t target, double theta) { push_2q_1p(c, PSI_GATE_CRX, control, target, theta); }
-void psi_apply_cry(struct PsiQuantumCircuit *c, size_t control, size_t target, double theta) { push_2q_1p(c, PSI_GATE_CRY, control, target, theta); }
-void psi_apply_crz(struct PsiQuantumCircuit *c, size_t control, size_t target, double theta) { push_2q_1p(c, PSI_GATE_CRZ, control, target, theta); }
-void psi_apply_cp(struct PsiQuantumCircuit *c, size_t control, size_t target, double theta) { push_2q_1p(c, PSI_GATE_CP, control, target, theta); }
+void psi_apply_cnot(struct PsiQuantumCircuit* c, size_t control, size_t target)
+{
+	push_2q(c, PSI_GATE_CNOT, control, target);
+}
+void psi_apply_cz(struct PsiQuantumCircuit* c, size_t control, size_t target)
+{
+	push_2q(c, PSI_GATE_CZ, control, target);
+}
+void psi_apply_swap(struct PsiQuantumCircuit* c, size_t qubit1, size_t qubit2)
+{
+	push_2q(c, PSI_GATE_SWAP, qubit1, qubit2);
+}
+void psi_apply_crx(struct PsiQuantumCircuit* c, size_t control, size_t target, double theta)
+{
+	push_2q_1p(c, PSI_GATE_CRX, control, target, theta);
+}
+void psi_apply_cry(struct PsiQuantumCircuit* c, size_t control, size_t target, double theta)
+{
+	push_2q_1p(c, PSI_GATE_CRY, control, target, theta);
+}
+void psi_apply_crz(struct PsiQuantumCircuit* c, size_t control, size_t target, double theta)
+{
+	push_2q_1p(c, PSI_GATE_CRZ, control, target, theta);
+}
+void psi_apply_cp(struct PsiQuantumCircuit* c, size_t control, size_t target, double theta)
+{
+	push_2q_1p(c, PSI_GATE_CP, control, target, theta);
+}
 
-void psi_apply_ccnot(struct PsiQuantumCircuit *c, size_t control1, size_t control2, size_t target) { push_3q(c, PSI_GATE_CCNOT, control1, control2, target); }
-void psi_apply_cswap(struct PsiQuantumCircuit *c, size_t control, size_t target1, size_t target2) { push_3q(c, PSI_GATE_CSWAP, control, target1, target2); }
+void psi_apply_ccnot(struct PsiQuantumCircuit* c, size_t control1, size_t control2, size_t target)
+{
+	push_3q(c, PSI_GATE_CCNOT, control1, control2, target);
+}
+void psi_apply_cswap(struct PsiQuantumCircuit* c, size_t control, size_t target1, size_t target2)
+{
+	push_3q(c, PSI_GATE_CSWAP, control, target1, target2);
+}
 
-void psi_measure(struct PsiQuantumCircuit *c, size_t qubit, size_t classical)
+void psi_measure(struct PsiQuantumCircuit* c, size_t qubit, size_t classical)
 {
 	if (classical >= c->num_classical)
 		c->num_classical = classical + 1;
 
-	struct PsiGateOp *op = append_op(c);
+	struct PsiGateOp* op = append_op(c);
 	op->kind = PSI_GATE_MEASURE;
 	op->qubits[0] = qubit;
 	op->qubit_count = 1;
 	op->classical = classical;
 }
 
-void psi_measure_all(struct PsiQuantumCircuit *c)
+void psi_measure_all(struct PsiQuantumCircuit* c)
 {
 	for (size_t i = 0; i < c->num_qubits; i++)
 		psi_measure(c, i, i);
 }
 
-void psi_apply_custom(struct PsiQuantumCircuit *c, struct PsiCustomGate gate, const size_t *targets, size_t count)
+void psi_apply_custom(struct PsiQuantumCircuit* c, struct PsiCustomGate gate, const size_t* targets,
+                      size_t count)
 {
-	struct PsiCustomGate *owned = malloc(sizeof(struct PsiCustomGate));
+	struct PsiCustomGate* owned = malloc(sizeof(struct PsiCustomGate));
 	assert(owned != NULL);
 	*owned = gate;
 
-	size_t *owned_targets = malloc(count * sizeof(size_t));
+	size_t* owned_targets = malloc(count * sizeof(size_t));
 	assert(owned_targets != NULL || count == 0);
 	if (count > 0)
 		memcpy(owned_targets, targets, count * sizeof(size_t));
 
-	struct PsiGateOp *op = append_op(c);
+	struct PsiGateOp* op = append_op(c);
 	op->kind = PSI_GATE_CUSTOM;
 	op->custom = owned;
 	op->custom_targets = owned_targets;
diff --git a/src/core/classical_components.c b/src/core/classical_components.c
index 836784c..3ffeb88 100644
--- a/src/core/classical_components.c
+++ b/src/core/classical_components.c
@@ -3,32 +3,31 @@
 #include <assert.h>
 #include <stdlib.h>
 
-struct PsiClassicalBit psi_new_classical_bit(const char *name, bool state)
+struct PsiClassicalBit psi_new_classical_bit(const char* name, bool state)
 {
-	return (struct PsiClassicalBit)
-	{
+	return (struct PsiClassicalBit){
 		name,
 		state,
 	};
 }
 
-struct PsiClassicalRegister psi_new_classical_register(const char *name, const char **names, size_t count)
+struct PsiClassicalRegister psi_new_classical_register(const char* name, const char** names,
+                                                       size_t count)
 {
-	struct PsiClassicalBit *bits = malloc(count * sizeof(struct PsiClassicalBit));
+	struct PsiClassicalBit* bits = malloc(count * sizeof(struct PsiClassicalBit));
 	assert(bits != NULL || count == 0);
 
 	for (size_t i = 0; i < count; i++)
 		bits[i] = psi_new_classical_bit(names[i], false);
 
-	return (struct PsiClassicalRegister)
-	{
+	return (struct PsiClassicalRegister){
 		name,
 		bits,
 		count,
 	};
 }
 
-void psi_free_classical_register(struct PsiClassicalRegister *reg)
+void psi_free_classical_register(struct PsiClassicalRegister* reg)
 {
 	free(reg->bits);
 	reg->bits = NULL;
diff --git a/src/core/custom_gate.c b/src/core/custom_gate.c
index 3b3a7e2..dd41635 100644
--- a/src/core/custom_gate.c
+++ b/src/core/custom_gate.c
@@ -7,7 +7,7 @@
 
 #include "core/gates.h"
 
-struct PsiCustomGate psi_new_custom_gate_from_matrix(const char *name, struct PsiMatrix matrix)
+struct PsiCustomGate psi_new_custom_gate_from_matrix(const char* name, struct PsiMatrix matrix)
 {
 	assert(matrix.rows == matrix.cols);
 
@@ -27,9 +27,11 @@ struct PsiCustomGate psi_new_custom_gate_from_matrix(const char *name, struct Ps
 	return gate;
 }
 
-struct PsiCustomGate psi_new_custom_gate_from_composite(const char *name, size_t num_qubits, const struct PsiCompositeGateOp *ops, size_t op_count)
+struct PsiCustomGate psi_new_custom_gate_from_composite(const char* name, size_t num_qubits,
+                                                        const struct PsiCompositeGateOp* ops,
+                                                        size_t op_count)
 {
-	struct PsiCompositeGateOp *owned = malloc(op_count * sizeof(struct PsiCompositeGateOp));
+	struct PsiCompositeGateOp* owned = malloc(op_count * sizeof(struct PsiCompositeGateOp));
 	assert(owned != NULL || op_count == 0);
 
 	if (op_count > 0)
@@ -45,7 +47,7 @@ struct PsiCustomGate psi_new_custom_gate_from_composite(const char *name, size_t
 	return gate;
 }
 
-void psi_free_custom_gate(struct PsiCustomGate *gate)
+void psi_free_custom_gate(struct PsiCustomGate* gate)
 {
 	if (gate->kind == PSI_CUSTOM_GATE_MATRIX)
 	{
@@ -78,7 +80,7 @@ static struct PsiQuantumGate op_gate(enum PsiCompositeOp op)
 	return psi_identity_gate();
 }
 
-static bool find_target(const size_t *targets, size_t count, size_t q, size_t *pos)
+static bool find_target(const size_t* targets, size_t count, size_t q, size_t* pos)
 {
 	for (size_t i = 0; i < count; i++)
 		if (targets[i] == q)
@@ -90,7 +92,8 @@ static bool find_target(const size_t *targets, size_t count, size_t q, size_t *p
 	return false;
 }
 
-static struct PsiMatrix build_full_operator(struct PsiMatrix gate_matrix, const size_t *targets, size_t num_gate_qubits, size_t total_qubits)
+static struct PsiMatrix build_full_operator(struct PsiMatrix gate_matrix, const size_t* targets,
+                                            size_t num_gate_qubits, size_t total_qubits)
 {
 	size_t dim = (size_t)1 << total_qubits;
 	size_t gate_dim = gate_matrix.rows;
@@ -139,7 +142,8 @@ static struct PsiMatrix compute_composite_matrix(struct PsiCustomGate gate)
 		struct PsiCompositeGateOp step = gate.definition.composite.ops[i];
 		struct PsiQuantumGate g = op_gate(step.op);
 
-		struct PsiMatrix full = build_full_operator(g.matrix, step.targets, step.target_count, gate.num_qubits);
+		struct PsiMatrix full =
+				build_full_operator(g.matrix, step.targets, step.target_count, gate.num_qubits);
 		struct PsiMatrix next = psi_dot_matrix(full, result);
 
 		psi_free_matrix(&full);
@@ -154,7 +158,8 @@ static struct PsiMatrix compute_composite_matrix(struct PsiCustomGate gate)
 struct PsiQuantumGate psi_to_quantum_gate(struct PsiCustomGate gate)
 {
 	if (gate.kind == PSI_CUSTOM_GATE_MATRIX)
-		return psi_new_quantum_gate(gate.name, psi_clone_matrix(gate.definition.matrix), gate.num_qubits);
+		return psi_new_quantum_gate(gate.name, psi_clone_matrix(gate.definition.matrix),
+		                            gate.num_qubits);
 
 	return psi_new_quantum_gate(gate.name, compute_composite_matrix(gate), gate.num_qubits);
 }
diff --git a/src/core/gates.c b/src/core/gates.c
index 0eb2997..51605c9 100644
--- a/src/core/gates.c
+++ b/src/core/gates.c
@@ -9,9 +9,8 @@ struct PsiMatrix psi_rx_matrix(double theta)
 	double cos_h = cos(theta / 2.0);
 	double sin_h = sin(theta / 2.0);
 
-	return psi_matrix(2, 2,
-		psi_new_complex(cos_h, 0.0), psi_new_complex(0.0, -sin_h),
-		psi_new_complex(0.0, -sin_h), psi_new_complex(cos_h, 0.0));
+	return psi_matrix(2, 2, psi_new_complex(cos_h, 0.0), psi_new_complex(0.0, -sin_h),
+	                  psi_new_complex(0.0, -sin_h), psi_new_complex(cos_h, 0.0));
 }
 
 struct PsiMatrix psi_ry_matrix(double theta)
@@ -19,25 +18,22 @@ struct PsiMatrix psi_ry_matrix(double theta)
 	double cos_h = cos(theta / 2.0);
 	double sin_h = sin(theta / 2.0);
 
-	return psi_matrix(2, 2,
-		psi_new_complex(cos_h, 0.0), psi_new_complex(-sin_h, 0.0),
-		psi_new_complex(sin_h, 0.0), psi_new_complex(cos_h, 0.0));
+	return psi_matrix(2, 2, psi_new_complex(cos_h, 0.0), psi_new_complex(-sin_h, 0.0),
+	                  psi_new_complex(sin_h, 0.0), psi_new_complex(cos_h, 0.0));
 }
 
 struct PsiMatrix psi_rz_matrix(double theta)
 {
 	double half = theta / 2.0;
 
-	return psi_matrix(2, 2,
-		psi_new_complex(cos(half), -sin(half)), psi_new_complex(0.0, 0.0),
-		psi_new_complex(0.0, 0.0), psi_new_complex(cos(half), sin(half)));
+	return psi_matrix(2, 2, psi_new_complex(cos(half), -sin(half)), psi_new_complex(0.0, 0.0),
+	                  psi_new_complex(0.0, 0.0), psi_new_complex(cos(half), sin(half)));
 }
 
 struct PsiMatrix psi_p_matrix(double theta)
 {
-	return psi_matrix(2, 2,
-		psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0),
-		psi_new_complex(0.0, 0.0), psi_new_complex(cos(theta), sin(theta)));
+	return psi_matrix(2, 2, psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0),
+	                  psi_new_complex(0.0, 0.0), psi_new_complex(cos(theta), sin(theta)));
 }
 
 struct PsiMatrix psi_u1_matrix(double lambda)
@@ -47,11 +43,11 @@ struct PsiMatrix psi_u1_matrix(double lambda)
 
 struct PsiMatrix psi_u2_matrix(double phi, double lambda)
 {
-	return psi_matrix(2, 2,
-		psi_new_complex(INV_SQRT_2, 0.0),
-		psi_new_complex(-INV_SQRT_2 * cos(lambda), -INV_SQRT_2 * sin(lambda)),
-		psi_new_complex(INV_SQRT_2 * cos(phi), INV_SQRT_2 * sin(phi)),
-		psi_new_complex(cos(phi + lambda) * INV_SQRT_2, sin(phi + lambda) * INV_SQRT_2));
+	return psi_matrix(
+			2, 2, psi_new_complex(INV_SQRT_2, 0.0),
+			psi_new_complex(-INV_SQRT_2 * cos(lambda), -INV_SQRT_2 * sin(lambda)),
+			psi_new_complex(INV_SQRT_2 * cos(phi), INV_SQRT_2 * sin(phi)),
+			psi_new_complex(cos(phi + lambda) * INV_SQRT_2, sin(phi + lambda) * INV_SQRT_2));
 }
 
 struct PsiMatrix psi_u3_matrix(double theta, double phi, double lambda)
@@ -59,11 +55,10 @@ struct PsiMatrix psi_u3_matrix(double theta, double phi, double lambda)
 	double cos_h = cos(theta / 2.0);
 	double sin_h = sin(theta / 2.0);
 
-	return psi_matrix(2, 2,
-		psi_new_complex(cos_h, 0.0),
-		psi_new_complex(-sin_h * cos(lambda), -sin_h * sin(lambda)),
-		psi_new_complex(sin_h * cos(phi), sin_h * sin(phi)),
-		psi_new_complex(cos_h * cos(phi + lambda), cos_h * sin(phi + lambda)));
+	return psi_matrix(2, 2, psi_new_complex(cos_h, 0.0),
+	                  psi_new_complex(-sin_h * cos(lambda), -sin_h * sin(lambda)),
+	                  psi_new_complex(sin_h * cos(phi), sin_h * sin(phi)),
+	                  psi_new_complex(cos_h * cos(phi + lambda), cos_h * sin(phi + lambda)));
 }
 
 struct PsiMatrix psi_crx_matrix(double theta)
@@ -115,90 +110,83 @@ struct PsiMatrix psi_cp_matrix(double theta)
 
 struct PsiQuantumGate psi_hadamard_gate(void)
 {
-	struct PsiMatrix m = psi_matrix(2, 2,
-		psi_new_complex(INV_SQRT_2, 0.0), psi_new_complex(INV_SQRT_2, 0.0),
-		psi_new_complex(INV_SQRT_2, 0.0), psi_new_complex(-INV_SQRT_2, 0.0));
+	struct PsiMatrix m =
+			psi_matrix(2, 2, psi_new_complex(INV_SQRT_2, 0.0), psi_new_complex(INV_SQRT_2, 0.0),
+			           psi_new_complex(INV_SQRT_2, 0.0), psi_new_complex(-INV_SQRT_2, 0.0));
 
 	return psi_new_quantum_gate_from_matrix("H", m);
 }
 
 struct PsiQuantumGate psi_pauli_x_gate(void)
 {
-	struct PsiMatrix m = 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 PsiMatrix m = 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));
 
 	return psi_new_quantum_gate_from_matrix("X", m);
 }
 
 struct PsiQuantumGate psi_pauli_y_gate(void)
 {
-	struct PsiMatrix m = psi_matrix(2, 2,
-		psi_new_complex(0.0, 0.0), psi_new_complex(0.0, -1.0),
-		psi_new_complex(0.0, 1.0), psi_new_complex(0.0, 0.0));
+	struct PsiMatrix m = psi_matrix(2, 2, psi_new_complex(0.0, 0.0), psi_new_complex(0.0, -1.0),
+	                                psi_new_complex(0.0, 1.0), psi_new_complex(0.0, 0.0));
 
 	return psi_new_quantum_gate_from_matrix("Y", m);
 }
 
 struct PsiQuantumGate psi_pauli_z_gate(void)
 {
-	struct PsiMatrix m = psi_matrix(2, 2,
-		psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0),
-		psi_new_complex(0.0, 0.0), psi_new_complex(-1.0, 0.0));
+	struct PsiMatrix m = psi_matrix(2, 2, psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0),
+	                                psi_new_complex(0.0, 0.0), psi_new_complex(-1.0, 0.0));
 
 	return psi_new_quantum_gate_from_matrix("Z", m);
 }
 
 struct PsiQuantumGate psi_s_gate(void)
 {
-	struct PsiMatrix m = psi_matrix(2, 2,
-		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, 1.0));
+	struct PsiMatrix m = psi_matrix(2, 2, 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, 1.0));
 
 	return psi_new_quantum_gate_from_matrix("S", m);
 }
 
 struct PsiQuantumGate psi_t_gate(void)
 {
-	struct PsiMatrix m = psi_matrix(2, 2,
-		psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0),
-		psi_new_complex(0.0, 0.0), psi_new_complex(INV_SQRT_2, INV_SQRT_2));
+	struct PsiMatrix m =
+			psi_matrix(2, 2, psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0),
+			           psi_new_complex(0.0, 0.0), psi_new_complex(INV_SQRT_2, INV_SQRT_2));
 
 	return psi_new_quantum_gate_from_matrix("T", m);
 }
 
 struct PsiQuantumGate psi_sdg_gate(void)
 {
-	struct PsiMatrix m = psi_matrix(2, 2,
-		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, -1.0));
+	struct PsiMatrix m = psi_matrix(2, 2, 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, -1.0));
 
 	return psi_new_quantum_gate_from_matrix("S†", m);
 }
 
 struct PsiQuantumGate psi_tdg_gate(void)
 {
-	struct PsiMatrix m = psi_matrix(2, 2,
-		psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0),
-		psi_new_complex(0.0, 0.0), psi_new_complex(INV_SQRT_2, -INV_SQRT_2));
+	struct PsiMatrix m =
+			psi_matrix(2, 2, psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0),
+			           psi_new_complex(0.0, 0.0), psi_new_complex(INV_SQRT_2, -INV_SQRT_2));
 
 	return psi_new_quantum_gate_from_matrix("T†", m);
 }
 
 struct PsiQuantumGate psi_sx_gate(void)
 {
-	struct PsiMatrix m = psi_matrix(2, 2,
-		psi_new_complex(0.5, 0.5), psi_new_complex(0.5, -0.5),
-		psi_new_complex(0.5, -0.5), psi_new_complex(0.5, 0.5));
+	struct PsiMatrix m = psi_matrix(2, 2, psi_new_complex(0.5, 0.5), psi_new_complex(0.5, -0.5),
+	                                psi_new_complex(0.5, -0.5), psi_new_complex(0.5, 0.5));
 
 	return psi_new_quantum_gate_from_matrix("√X", m);
 }
 
 struct PsiQuantumGate psi_sxdg_gate(void)
 {
-	struct PsiMatrix m = psi_matrix(2, 2,
-		psi_new_complex(0.5, -0.5), psi_new_complex(0.5, 0.5),
-		psi_new_complex(0.5, 0.5), psi_new_complex(0.5, -0.5));
+	struct PsiMatrix m = psi_matrix(2, 2, psi_new_complex(0.5, -0.5), psi_new_complex(0.5, 0.5),
+	                                psi_new_complex(0.5, 0.5), psi_new_complex(0.5, -0.5));
 
 	return psi_new_quantum_gate_from_matrix("√X†", m);
 }
diff --git a/src/core/kernel.c b/src/core/kernel.c
index f8f7dfd..def0782 100644
--- a/src/core/kernel.c
+++ b/src/core/kernel.c
@@ -6,19 +6,19 @@
 #include <stdlib.h>
 #include <string.h>
 
-static char *dup_string(const char *s)
+static char* dup_string(const char* s)
 {
 	size_t n = strlen(s) + 1;
-	char *p = malloc(n);
+	char* p = malloc(n);
 	assert(p != NULL);
 	memcpy(p, s, n);
 
 	return p;
 }
 
-static size_t *dup_targets(const size_t *targets, size_t count)
+static size_t* dup_targets(const size_t* targets, size_t count)
 {
-	size_t *p = malloc(count * sizeof(size_t));
+	size_t* p = malloc(count * sizeof(size_t));
 	assert(p != NULL || count == 0);
 
 	if (count > 0)
@@ -27,29 +27,30 @@ static size_t *dup_targets(const size_t *targets, size_t count)
 	return p;
 }
 
-static bool starts_with(const char *s, const char *prefix)
+static bool starts_with(const char* s, const char* prefix)
 {
 	return strncmp(s, prefix, strlen(prefix)) == 0;
 }
 
-static enum PsiGateType detect_gate_type(const char *name, struct PsiMatrix matrix)
+static enum PsiGateType detect_gate_type(const char* name, struct PsiMatrix matrix)
 {
-	static const char *diagonal[] = { "Z", "S", "T", "Sdg", "Tdg", "Rz", "P", "U1", "CZ", "CP", "CRz" };
+	static const char* diagonal[] = { "Z", "S",  "T",  "Sdg", "Tdg", "Rz",
+	                                  "P", "U1", "CZ", "CP",  "CRz" };
 	for (size_t i = 0; i < sizeof(diagonal) / sizeof(diagonal[0]); i++)
 		if (starts_with(name, diagonal[i]))
 			return PSI_GATE_TYPE_DIAGONAL;
 
-	static const char *controlled[] = { "CNOT", "CZ", "SWAP", "CRx", "CRy", "CRz", "CP", "CCNOT", "CSWAP" };
+	static const char* controlled[] = { "CNOT", "CZ", "SWAP",  "CRx",  "CRy",
+	                                    "CRz",  "CP", "CCNOT", "CSWAP" };
 	for (size_t i = 0; i < sizeof(controlled) / sizeof(controlled[0]); i++)
 		if (starts_with(name, controlled[i]))
 			return PSI_GATE_TYPE_CONTROLLED;
 
 	if (matrix.rows == 2 && matrix.cols == 2)
 	{
-		bool is_diag = fabs(matrix.data[1].real) < 1e-10
-			&& fabs(matrix.data[1].imaginary) < 1e-10
-			&& fabs(matrix.data[2].real) < 1e-10
-			&& fabs(matrix.data[2].imaginary) < 1e-10;
+		bool is_diag = fabs(matrix.data[1].real) < 1e-10 &&
+				fabs(matrix.data[1].imaginary) < 1e-10 && fabs(matrix.data[2].real) < 1e-10 &&
+				fabs(matrix.data[2].imaginary) < 1e-10;
 		if (is_diag)
 			return PSI_GATE_TYPE_DIAGONAL;
 	}
@@ -57,7 +58,8 @@ static enum PsiGateType detect_gate_type(const char *name, struct PsiMatrix matr
 	return PSI_GATE_TYPE_NON_DIAGONAL;
 }
 
-struct PsiKernel psi_new_kernel(const char *name, struct PsiMatrix matrix, const size_t *targets, size_t target_count)
+struct PsiKernel psi_new_kernel(const char* name, struct PsiMatrix matrix, const size_t* targets,
+                                size_t target_count)
 {
 	struct PsiKernel kernel;
 	kernel.matrix = matrix;
@@ -81,7 +83,7 @@ struct PsiKernel psi_clone_kernel(struct PsiKernel kernel)
 	return copy;
 }
 
-void psi_free_kernel(struct PsiKernel *kernel)
+void psi_free_kernel(struct PsiKernel* kernel)
 {
 	psi_free_matrix(&kernel->matrix);
 	free(kernel->targets);
@@ -118,7 +120,8 @@ bool psi_kernels_commute(struct PsiKernel a, struct PsiKernel b)
 	if (!psi_kernels_share_qubits(a, b))
 		return true;
 
-	if (a.gate_type == PSI_GATE_TYPE_DIAGONAL && b.gate_type == PSI_GATE_TYPE_DIAGONAL && targets_equal(a, b))
+	if (a.gate_type == PSI_GATE_TYPE_DIAGONAL && b.gate_type == PSI_GATE_TYPE_DIAGONAL &&
+	    targets_equal(a, b))
 		return true;
 
 	return false;
@@ -132,17 +135,18 @@ bool psi_kernels_can_fuse(struct PsiKernel a, struct PsiKernel b)
 	return a.targets[0] == b.targets[0];
 }
 
-bool psi_fuse_kernels(struct PsiKernel a, struct PsiKernel b, struct PsiKernel *out)
+bool psi_fuse_kernels(struct PsiKernel a, struct PsiKernel b, struct PsiKernel* out)
 {
 	if (!psi_kernels_can_fuse(a, b))
 		return false;
 
-	enum PsiGateType new_type = a.gate_type == PSI_GATE_TYPE_DIAGONAL && b.gate_type == PSI_GATE_TYPE_DIAGONAL
-		? PSI_GATE_TYPE_DIAGONAL
-		: PSI_GATE_TYPE_NON_DIAGONAL;
+	enum PsiGateType new_type =
+			a.gate_type == PSI_GATE_TYPE_DIAGONAL && b.gate_type == PSI_GATE_TYPE_DIAGONAL
+			? PSI_GATE_TYPE_DIAGONAL
+			: PSI_GATE_TYPE_NON_DIAGONAL;
 
 	size_t name_len = strlen(a.name) + strlen(b.name) + 2;
-	char *fused_name = malloc(name_len);
+	char* fused_name = malloc(name_len);
 	assert(fused_name != NULL);
 	snprintf(fused_name, name_len, "%s+%s", a.name, b.name);
 
@@ -155,13 +159,14 @@ bool psi_fuse_kernels(struct PsiKernel a, struct PsiKernel b, struct PsiKernel *
 	return true;
 }
 
-static struct PsiComplex *apply_kernel(const struct PsiComplex *state, struct PsiKernel kernel, size_t num_qubits)
+static struct PsiComplex* apply_kernel(const struct PsiComplex* state, struct PsiKernel kernel,
+                                       size_t num_qubits)
 {
 	size_t dim = (size_t)1 << num_qubits;
 	size_t g = kernel.target_count;
 	size_t gate_dim = (size_t)1 << g;
 
-	size_t *target_bits = malloc(g * sizeof(size_t));
+	size_t* target_bits = malloc(g * sizeof(size_t));
 	assert(target_bits != NULL || g == 0);
 	for (size_t k = 0; k < g; k++)
 		target_bits[k] = num_qubits - 1 - kernel.targets[k];
@@ -170,7 +175,7 @@ static struct PsiComplex *apply_kernel(const struct PsiComplex *state, struct Ps
 	for (size_t k = 0; k < g; k++)
 		non_target_mask &= ~((size_t)1 << target_bits[k]);
 
-	struct PsiComplex *new_state = malloc(dim * sizeof(struct PsiComplex));
+	struct PsiComplex* new_state = malloc(dim * sizeof(struct PsiComplex));
 	assert(new_state != NULL);
 
 	for (size_t i = 0; i < dim; i++)
@@ -213,7 +218,7 @@ struct PsiKernelBatch psi_new_kernel_batch(size_t num_qubits)
 	return batch;
 }
 
-void psi_free_kernel_batch(struct PsiKernelBatch *batch)
+void psi_free_kernel_batch(struct PsiKernelBatch* batch)
 {
 	for (size_t i = 0; i < batch->count; i++)
 		psi_free_kernel(&batch->kernels[i]);
@@ -224,7 +229,7 @@ void psi_free_kernel_batch(struct PsiKernelBatch *batch)
 	batch->capacity = 0;
 }
 
-void psi_add_kernel(struct PsiKernelBatch *batch, struct PsiKernel kernel)
+void psi_add_kernel(struct PsiKernelBatch* batch, struct PsiKernel kernel)
 {
 	if (batch->count == batch->capacity)
 	{
@@ -237,13 +242,13 @@ void psi_add_kernel(struct PsiKernelBatch *batch, struct PsiKernel kernel)
 	batch->kernels[batch->count++] = kernel;
 }
 
-void psi_optimize_kernel_batch(struct PsiKernelBatch *batch)
+void psi_optimize_kernel_batch(struct PsiKernelBatch* batch)
 {
 	if (batch->count < 2)
 		return;
 
 	size_t original = batch->count;
-	struct PsiKernel *out = malloc(original * sizeof(struct PsiKernel));
+	struct PsiKernel* out = malloc(original * sizeof(struct PsiKernel));
 	assert(out != NULL);
 	size_t out_count = 0;
 
@@ -273,27 +278,28 @@ void psi_optimize_kernel_batch(struct PsiKernelBatch *batch)
 	batch->capacity = original;
 }
 
-void psi_execute_kernel_batch(struct PsiKernelBatch batch, struct PsiVector *state)
+void psi_execute_kernel_batch(struct PsiKernelBatch batch, struct PsiVector* state)
 {
 	size_t dim = (size_t)1 << batch.num_qubits;
 	assert(state->size == dim);
 
 	for (size_t i = 0; i < batch.count; i++)
 	{
-		struct PsiComplex *next = apply_kernel(state->data, batch.kernels[i], batch.num_qubits);
+		struct PsiComplex* next = apply_kernel(state->data, batch.kernels[i], batch.num_qubits);
 		free(state->data);
 		state->data = next;
 	}
 }
 
-void psi_apply_kernel(struct PsiVector *state, struct PsiKernel kernel, size_t num_qubits)
+void psi_apply_kernel(struct PsiVector* state, struct PsiKernel kernel, size_t num_qubits)
 {
-	struct PsiComplex *next = apply_kernel(state->data, kernel, num_qubits);
+	struct PsiComplex* next = apply_kernel(state->data, kernel, num_qubits);
 	free(state->data);
 	state->data = next;
 }
 
-static void push_kernel(struct PsiKernel **kernels, size_t *count, size_t *capacity, struct PsiKernel kernel)
+static void push_kernel(struct PsiKernel** kernels, size_t* count, size_t* capacity,
+                        struct PsiKernel kernel)
 {
 	if (*count == *capacity)
 	{
@@ -315,7 +321,7 @@ static bool layer_can_add(struct PsiExecutionLayer layer, struct PsiKernel kerne
 	return true;
 }
 
-static void free_layer(struct PsiExecutionLayer *layer)
+static void free_layer(struct PsiExecutionLayer* layer)
 {
 	for (size_t i = 0; i < layer->count; i++)
 		psi_free_kernel(&layer->kernels[i]);
@@ -341,7 +347,7 @@ struct PsiStructureAwareBatch psi_new_structure_aware_batch(size_t num_qubits)
 	return batch;
 }
 
-static void clear_layers(struct PsiStructureAwareBatch *batch)
+static void clear_layers(struct PsiStructureAwareBatch* batch)
 {
 	for (size_t i = 0; i < batch->layer_count; i++)
 		free_layer(&batch->layers[i]);
@@ -352,7 +358,7 @@ static void clear_layers(struct PsiStructureAwareBatch *batch)
 	batch->layer_capacity = 0;
 }
 
-void psi_free_structure_aware_batch(struct PsiStructureAwareBatch *batch)
+void psi_free_structure_aware_batch(struct PsiStructureAwareBatch* batch)
 {
 	for (size_t i = 0; i < batch->count; i++)
 		psi_free_kernel(&batch->kernels[i]);
@@ -364,13 +370,13 @@ void psi_free_structure_aware_batch(struct PsiStructureAwareBatch *batch)
 	clear_layers(batch);
 }
 
-void psi_add_structure_aware_kernel(struct PsiStructureAwareBatch *batch, struct PsiKernel kernel)
+void psi_add_structure_aware_kernel(struct PsiStructureAwareBatch* batch, struct PsiKernel kernel)
 {
 	push_kernel(&batch->kernels, &batch->count, &batch->capacity, kernel);
 	batch->optimised = false;
 }
 
-static struct PsiKernel remove_kernel_at(struct PsiStructureAwareBatch *batch, size_t index)
+static struct PsiKernel remove_kernel_at(struct PsiStructureAwareBatch* batch, size_t index)
 {
 	struct PsiKernel removed = batch->kernels[index];
 	for (size_t i = index; i + 1 < batch->count; i++)
@@ -380,7 +386,8 @@ static struct PsiKernel remove_kernel_at(struct PsiStructureAwareBatch *batch, s
 	return removed;
 }
 
-static void insert_kernel_at(struct PsiStructureAwareBatch *batch, size_t index, struct PsiKernel kernel)
+static void insert_kernel_at(struct PsiStructureAwareBatch* batch, size_t index,
+                             struct PsiKernel kernel)
 {
 	if (batch->count == batch->capacity)
 	{
@@ -397,7 +404,7 @@ static void insert_kernel_at(struct PsiStructureAwareBatch *batch, size_t index,
 	batch->count++;
 }
 
-static void reorder_commuting_gates(struct PsiStructureAwareBatch *batch)
+static void reorder_commuting_gates(struct PsiStructureAwareBatch* batch)
 {
 	bool changed = true;
 	size_t iterations = 0;
@@ -413,9 +420,8 @@ static void reorder_commuting_gates(struct PsiStructureAwareBatch *batch)
 			struct PsiKernel current = batch->kernels[i];
 			struct PsiKernel next = batch->kernels[i + 1];
 
-			if (current.target_count != 1 || next.target_count != 1
-				|| current.targets[0] == next.targets[0]
-				|| !psi_kernels_commute(current, next))
+			if (current.target_count != 1 || next.target_count != 1 ||
+			    current.targets[0] == next.targets[0] || !psi_kernels_commute(current, next))
 				continue;
 
 			for (size_t j = i + 2; j < batch->count; j++)
@@ -428,7 +434,8 @@ static void reorder_commuting_gates(struct PsiStructureAwareBatch *batch)
 				for (size_t k = i + 1; k < j; k++)
 				{
 					struct PsiKernel between = batch->kernels[k];
-					if (psi_kernels_share_qubits(between, current) && !psi_kernels_commute(current, between))
+					if (psi_kernels_share_qubits(between, current) &&
+					    !psi_kernels_commute(current, between))
 					{
 						can_move = false;
 						break;
@@ -447,7 +454,7 @@ static void reorder_commuting_gates(struct PsiStructureAwareBatch *batch)
 	}
 }
 
-static void multi_pass_fusion(struct PsiStructureAwareBatch *batch)
+static void multi_pass_fusion(struct PsiStructureAwareBatch* batch)
 {
 	bool changed = true;
 	size_t iterations = 0;
@@ -458,7 +465,7 @@ static void multi_pass_fusion(struct PsiStructureAwareBatch *batch)
 		changed = false;
 		iterations++;
 
-		struct PsiKernel *new_kernels = NULL;
+		struct PsiKernel* new_kernels = NULL;
 		size_t new_count = 0;
 		size_t new_capacity = 0;
 
@@ -490,7 +497,7 @@ static void multi_pass_fusion(struct PsiStructureAwareBatch *batch)
 	}
 }
 
-static void build_execution_layers(struct PsiStructureAwareBatch *batch)
+static void build_execution_layers(struct PsiStructureAwareBatch* batch)
 {
 	clear_layers(batch);
 
@@ -502,8 +509,9 @@ static void build_execution_layers(struct PsiStructureAwareBatch *batch)
 		for (size_t l = 0; l < batch->layer_count; l++)
 			if (layer_can_add(batch->layers[l], kernel))
 			{
-				struct PsiExecutionLayer *layer = &batch->layers[l];
-				push_kernel(&layer->kernels, &layer->count, &layer->capacity, psi_clone_kernel(kernel));
+				struct PsiExecutionLayer* layer = &batch->layers[l];
+				push_kernel(&layer->kernels, &layer->count, &layer->capacity,
+				            psi_clone_kernel(kernel));
 				placed = true;
 				break;
 			}
@@ -528,7 +536,7 @@ static void build_execution_layers(struct PsiStructureAwareBatch *batch)
 	}
 }
 
-void psi_optimize_structure_aware_batch(struct PsiStructureAwareBatch *batch)
+void psi_optimize_structure_aware_batch(struct PsiStructureAwareBatch* batch)
 {
 	if (batch->optimised || batch->count < 2)
 		return;
@@ -539,20 +547,21 @@ void psi_optimize_structure_aware_batch(struct PsiStructureAwareBatch *batch)
 	batch->optimised = true;
 }
 
-void psi_execute_structure_aware_batch(struct PsiStructureAwareBatch batch, struct PsiVector *state)
+void psi_execute_structure_aware_batch(struct PsiStructureAwareBatch batch, struct PsiVector* state)
 {
 	size_t dim = (size_t)1 << batch.num_qubits;
 	assert(state->size == dim);
 
 	for (size_t i = 0; i < batch.count; i++)
 	{
-		struct PsiComplex *next = apply_kernel(state->data, batch.kernels[i], batch.num_qubits);
+		struct PsiComplex* next = apply_kernel(state->data, batch.kernels[i], batch.num_qubits);
 		free(state->data);
 		state->data = next;
 	}
 }
 
-void psi_execute_structure_aware_batch_layered(struct PsiStructureAwareBatch batch, struct PsiVector *state)
+void psi_execute_structure_aware_batch_layered(struct PsiStructureAwareBatch batch,
+                                               struct PsiVector* state)
 {
 	size_t dim = (size_t)1 << batch.num_qubits;
 	assert(state->size == dim);
@@ -560,7 +569,8 @@ void psi_execute_structure_aware_batch_layered(struct PsiStructureAwareBatch bat
 	for (size_t l = 0; l < batch.layer_count; l++)
 		for (size_t k = 0; k < batch.layers[l].count; k++)
 		{
-			struct PsiComplex *next = apply_kernel(state->data, batch.layers[l].kernels[k], batch.num_qubits);
+			struct PsiComplex* next =
+					apply_kernel(state->data, batch.layers[l].kernels[k], batch.num_qubits);
 			free(state->data);
 			state->data = next;
 		}
diff --git a/src/core/noise.c b/src/core/noise.c
index 9aaf1d9..5f141c8 100644
--- a/src/core/noise.c
+++ b/src/core/noise.c
@@ -5,30 +5,30 @@
 #include <stdlib.h>
 #include <string.h>
 
-struct PsiKrausOperator psi_new_kraus_operator(const char *name, struct PsiMatrix matrix)
+struct PsiKrausOperator psi_new_kraus_operator(const char* name, struct PsiMatrix matrix)
 {
-	return (struct PsiKrausOperator)
-	{
+	return (struct PsiKrausOperator){
 		matrix,
 		name,
 	};
 }
 
-void psi_free_kraus_operator(struct PsiKrausOperator *op)
+void psi_free_kraus_operator(struct PsiKrausOperator* op)
 {
 	psi_free_matrix(&op->matrix);
 }
 
-struct PsiNoiseChannel psi_new_noise_channel(const char *name, const struct PsiKrausOperator *operators, size_t count, size_t num_qubits)
+struct PsiNoiseChannel psi_new_noise_channel(const char* name,
+                                             const struct PsiKrausOperator* operators, size_t count,
+                                             size_t num_qubits)
 {
-	struct PsiKrausOperator *owned = malloc(count * sizeof(struct PsiKrausOperator));
+	struct PsiKrausOperator* owned = malloc(count * sizeof(struct PsiKrausOperator));
 	assert(owned != NULL || count == 0);
 
 	if (count > 0)
 		memcpy(owned, operators, count * sizeof(struct PsiKrausOperator));
 
-	return (struct PsiNoiseChannel)
-	{
+	return (struct PsiNoiseChannel){
 		name,
 		owned,
 		count,
@@ -36,7 +36,7 @@ struct PsiNoiseChannel psi_new_noise_channel(const char *name, const struct PsiK
 	};
 }
 
-void psi_free_noise_channel(struct PsiNoiseChannel *channel)
+void psi_free_noise_channel(struct PsiNoiseChannel* channel)
 {
 	for (size_t i = 0; i < channel->operator_count; i++)
 		psi_free_matrix(&channel->operators[i].matrix);
@@ -52,18 +52,22 @@ struct PsiNoiseChannel psi_depolarising_channel(double p)
 	double sqrt_p3 = sqrt(p / 3.0);
 
 	struct PsiKrausOperator ops[] = {
-		psi_new_kraus_operator("K0", psi_matrix(2, 2,
-			psi_new_complex(sqrt_1_p, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_1_p, 0.0))),
-		psi_new_kraus_operator("K1(X)", psi_matrix(2, 2,
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_p3, 0.0),
-			psi_new_complex(sqrt_p3, 0.0), psi_new_complex(0.0, 0.0))),
-		psi_new_kraus_operator("K2(Y)", psi_matrix(2, 2,
-			psi_new_complex(0.0, 0.0), psi_new_complex(0.0, -sqrt_p3),
-			psi_new_complex(0.0, sqrt_p3), psi_new_complex(0.0, 0.0))),
-		psi_new_kraus_operator("K3(Z)", psi_matrix(2, 2,
-			psi_new_complex(sqrt_p3, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(-sqrt_p3, 0.0))),
+		psi_new_kraus_operator("K0",
+		                       psi_matrix(2, 2, psi_new_complex(sqrt_1_p, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_1_p, 0.0))),
+		psi_new_kraus_operator(
+				"K1(X)",
+				psi_matrix(2, 2, psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_p3, 0.0),
+				           psi_new_complex(sqrt_p3, 0.0), psi_new_complex(0.0, 0.0))),
+		psi_new_kraus_operator(
+				"K2(Y)",
+				psi_matrix(2, 2, psi_new_complex(0.0, 0.0), psi_new_complex(0.0, -sqrt_p3),
+				           psi_new_complex(0.0, sqrt_p3), psi_new_complex(0.0, 0.0))),
+		psi_new_kraus_operator("K3(Z)",
+		                       psi_matrix(2, 2, psi_new_complex(sqrt_p3, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(-sqrt_p3, 0.0))),
 	};
 
 	return psi_new_noise_channel("Depolarising", ops, 4, 1);
@@ -75,12 +79,14 @@ struct PsiNoiseChannel psi_amplitude_damping_channel(double gamma)
 	double sqrt_1_gamma = sqrt(1.0 - gamma);
 
 	struct PsiKrausOperator ops[] = {
-		psi_new_kraus_operator("K0", psi_matrix(2, 2,
-			psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_1_gamma, 0.0))),
-		psi_new_kraus_operator("K1", psi_matrix(2, 2,
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_gamma, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0))),
+		psi_new_kraus_operator("K0",
+		                       psi_matrix(2, 2, psi_new_complex(1.0, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_1_gamma, 0.0))),
+		psi_new_kraus_operator("K1",
+		                       psi_matrix(2, 2, psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_gamma, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0))),
 	};
 
 	return psi_new_noise_channel("AmplitudeDamping", ops, 2, 1);
@@ -92,12 +98,14 @@ struct PsiNoiseChannel psi_phase_damping_channel(double gamma)
 	double sqrt_1_gamma = sqrt(1.0 - gamma);
 
 	struct PsiKrausOperator ops[] = {
-		psi_new_kraus_operator("K0", psi_matrix(2, 2,
-			psi_new_complex(1.0, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_1_gamma, 0.0))),
-		psi_new_kraus_operator("K1", psi_matrix(2, 2,
-			psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_gamma, 0.0))),
+		psi_new_kraus_operator("K0",
+		                       psi_matrix(2, 2, psi_new_complex(1.0, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_1_gamma, 0.0))),
+		psi_new_kraus_operator("K1",
+		                       psi_matrix(2, 2, psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_gamma, 0.0))),
 	};
 
 	return psi_new_noise_channel("PhaseDamping", ops, 2, 1);
@@ -109,12 +117,14 @@ struct PsiNoiseChannel psi_bit_flip_channel(double p)
 	double sqrt_p = sqrt(p);
 
 	struct PsiKrausOperator ops[] = {
-		psi_new_kraus_operator("K0(I)", psi_matrix(2, 2,
-			psi_new_complex(sqrt_1_p, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_1_p, 0.0))),
-		psi_new_kraus_operator("K1(X)", psi_matrix(2, 2,
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_p, 0.0),
-			psi_new_complex(sqrt_p, 0.0), psi_new_complex(0.0, 0.0))),
+		psi_new_kraus_operator("K0(I)",
+		                       psi_matrix(2, 2, psi_new_complex(sqrt_1_p, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_1_p, 0.0))),
+		psi_new_kraus_operator("K1(X)",
+		                       psi_matrix(2, 2, psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_p, 0.0),
+		                                  psi_new_complex(sqrt_p, 0.0), psi_new_complex(0.0, 0.0))),
 	};
 
 	return psi_new_noise_channel("BitFlip", ops, 2, 1);
@@ -126,12 +136,14 @@ struct PsiNoiseChannel psi_phase_flip_channel(double p)
 	double sqrt_p = sqrt(p);
 
 	struct PsiKrausOperator ops[] = {
-		psi_new_kraus_operator("K0(I)", psi_matrix(2, 2,
-			psi_new_complex(sqrt_1_p, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_1_p, 0.0))),
-		psi_new_kraus_operator("K1(Z)", psi_matrix(2, 2,
-			psi_new_complex(sqrt_p, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(-sqrt_p, 0.0))),
+		psi_new_kraus_operator("K0(I)",
+		                       psi_matrix(2, 2, psi_new_complex(sqrt_1_p, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_1_p, 0.0))),
+		psi_new_kraus_operator("K1(Z)",
+		                       psi_matrix(2, 2, psi_new_complex(sqrt_p, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(-sqrt_p, 0.0))),
 	};
 
 	return psi_new_noise_channel("PhaseFlip", ops, 2, 1);
@@ -143,12 +155,14 @@ struct PsiNoiseChannel psi_bit_phase_flip_channel(double p)
 	double sqrt_p = sqrt(p);
 
 	struct PsiKrausOperator ops[] = {
-		psi_new_kraus_operator("K0(I)", psi_matrix(2, 2,
-			psi_new_complex(sqrt_1_p, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_1_p, 0.0))),
-		psi_new_kraus_operator("K1(Y)", psi_matrix(2, 2,
-			psi_new_complex(0.0, 0.0), psi_new_complex(0.0, -sqrt_p),
-			psi_new_complex(0.0, sqrt_p), psi_new_complex(0.0, 0.0))),
+		psi_new_kraus_operator("K0(I)",
+		                       psi_matrix(2, 2, psi_new_complex(sqrt_1_p, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_1_p, 0.0))),
+		psi_new_kraus_operator("K1(Y)",
+		                       psi_matrix(2, 2, psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(0.0, -sqrt_p),
+		                                  psi_new_complex(0.0, sqrt_p), psi_new_complex(0.0, 0.0))),
 	};
 
 	return psi_new_noise_channel("BitPhaseFlip", ops, 2, 1);
@@ -162,18 +176,22 @@ struct PsiNoiseChannel psi_generalised_amplitude_damping_channel(double p, doubl
 	double sqrt_1_gamma = sqrt(1.0 - gamma);
 
 	struct PsiKrausOperator ops[] = {
-		psi_new_kraus_operator("K0", psi_matrix(2, 2,
-			psi_new_complex(sqrt_p, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_p * sqrt_1_gamma, 0.0))),
-		psi_new_kraus_operator("K1", psi_matrix(2, 2,
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_p * sqrt_gamma, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0))),
-		psi_new_kraus_operator("K2", psi_matrix(2, 2,
-			psi_new_complex(sqrt_1_p * sqrt_1_gamma, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(0.0, 0.0), psi_new_complex(sqrt_1_p, 0.0))),
-		psi_new_kraus_operator("K3", psi_matrix(2, 2,
-			psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
-			psi_new_complex(sqrt_1_p * sqrt_gamma, 0.0), psi_new_complex(0.0, 0.0))),
+		psi_new_kraus_operator("K0",
+		                       psi_matrix(2, 2, psi_new_complex(sqrt_p, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_p * sqrt_1_gamma, 0.0))),
+		psi_new_kraus_operator("K1",
+		                       psi_matrix(2, 2, psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_p * sqrt_gamma, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0))),
+		psi_new_kraus_operator("K2",
+		                       psi_matrix(2, 2, psi_new_complex(sqrt_1_p * sqrt_1_gamma, 0.0),
+		                                  psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+		                                  psi_new_complex(sqrt_1_p, 0.0))),
+		psi_new_kraus_operator(
+				"K3",
+				psi_matrix(2, 2, psi_new_complex(0.0, 0.0), psi_new_complex(0.0, 0.0),
+				           psi_new_complex(sqrt_1_p * sqrt_gamma, 0.0), psi_new_complex(0.0, 0.0))),
 	};
 
 	return psi_new_noise_channel("GeneralisedAmplitudeDamping", ops, 4, 1);
@@ -182,41 +200,40 @@ struct PsiNoiseChannel psi_generalised_amplitude_damping_channel(double p, doubl
 struct PsiDensityMatrix psi_new_density_matrix(size_t num_qubits)
 {
 	size_t dim = (size_t)1 << num_qubits;
-	struct PsiComplex *data = calloc(dim * dim, sizeof(struct PsiComplex));
+	struct PsiComplex* data = calloc(dim * dim, sizeof(struct PsiComplex));
 	assert(data != NULL);
 	data[0] = psi_new_complex(1.0, 0.0);
 
-	return (struct PsiDensityMatrix)
-	{
+	return (struct PsiDensityMatrix){
 		data,
 		dim,
 		num_qubits,
 	};
 }
 
-struct PsiDensityMatrix psi_new_density_matrix_from_state(const struct PsiComplex *state, size_t len)
+struct PsiDensityMatrix psi_new_density_matrix_from_state(const struct PsiComplex* state,
+                                                          size_t len)
 {
 	size_t dim = len;
 	size_t num_qubits = 0;
 	while (((size_t)1 << num_qubits) < dim)
 		num_qubits++;
 
-	struct PsiComplex *data = malloc(dim * dim * sizeof(struct PsiComplex));
+	struct PsiComplex* data = malloc(dim * dim * sizeof(struct PsiComplex));
 	assert(data != NULL);
 
 	for (size_t i = 0; i < dim; i++)
 		for (size_t j = 0; j < dim; j++)
 			data[i * dim + j] = psi_mul_complex(state[i], psi_conjugate_complex(state[j]));
 
-	return (struct PsiDensityMatrix)
-	{
+	return (struct PsiDensityMatrix){
 		data,
 		dim,
 		num_qubits,
 	};
 }
 
-void psi_free_density_matrix(struct PsiDensityMatrix *dm)
+void psi_free_density_matrix(struct PsiDensityMatrix* dm)
 {
 	free(dm->data);
 	dm->data = NULL;
@@ -230,7 +247,8 @@ struct PsiComplex psi_get_density_matrix(struct PsiDensityMatrix dm, size_t row,
 	return dm.data[row * dm.dim + col];
 }
 
-void psi_set_density_matrix(struct PsiDensityMatrix *dm, size_t row, size_t col, struct PsiComplex value)
+void psi_set_density_matrix(struct PsiDensityMatrix* dm, size_t row, size_t col,
+                            struct PsiComplex value)
 {
 	assert(row < dm->dim && col < dm->dim);
 	dm->data[row * dm->dim + col] = value;
@@ -250,7 +268,8 @@ double psi_purity_density_matrix(struct PsiDensityMatrix dm)
 	struct PsiComplex sum = psi_new_complex(0.0, 0.0);
 	for (size_t i = 0; i < dm.dim; i++)
 		for (size_t j = 0; j < dm.dim; j++)
-			sum = psi_add_complex(sum, psi_mul_complex(dm.data[i * dm.dim + j], dm.data[j * dm.dim + i]));
+			sum = psi_add_complex(
+					sum, psi_mul_complex(dm.data[i * dm.dim + j], dm.data[j * dm.dim + i]));
 
 	return sum.real;
 }
@@ -260,19 +279,20 @@ bool psi_is_pure_density_matrix(struct PsiDensityMatrix dm, double tolerance)
 	return fabs(psi_purity_density_matrix(dm) - 1.0) < tolerance;
 }
 
-void psi_density_matrix_probabilities(struct PsiDensityMatrix dm, double *out)
+void psi_density_matrix_probabilities(struct PsiDensityMatrix dm, double* out)
 {
 	for (size_t i = 0; i < dm.dim; i++)
 		out[i] = dm.data[i * dm.dim + i].real;
 }
 
-void psi_apply_unitary_density_matrix(struct PsiDensityMatrix *dm, struct PsiMatrix gate, const size_t *targets, size_t target_count)
+void psi_apply_unitary_density_matrix(struct PsiDensityMatrix* dm, struct PsiMatrix gate,
+                                      const size_t* targets, size_t target_count)
 {
 	size_t g = target_count;
 	size_t gate_dim = (size_t)1 << g;
 	size_t dim = dm->dim;
 
-	size_t *target_bits = malloc(g * sizeof(size_t));
+	size_t* target_bits = malloc(g * sizeof(size_t));
 	assert(target_bits != NULL || g == 0);
 	for (size_t t = 0; t < g; t++)
 		target_bits[t] = dm->num_qubits - 1 - targets[t];
@@ -281,7 +301,7 @@ void psi_apply_unitary_density_matrix(struct PsiDensityMatrix *dm, struct PsiMat
 	for (size_t t = 0; t < g; t++)
 		non_target_mask &= ~((size_t)1 << target_bits[t]);
 
-	struct PsiComplex *new_data = calloc(dim * dim, sizeof(struct PsiComplex));
+	struct PsiComplex* new_data = calloc(dim * dim, sizeof(struct PsiComplex));
 	assert(new_data != NULL);
 
 	for (size_t i = 0; i < dim; i++)
@@ -314,10 +334,12 @@ void psi_apply_unitary_density_matrix(struct PsiDensityMatrix *dm, struct PsiMat
 					}
 
 					struct PsiComplex u_ik = gate.data[tgt_i * gate_dim + k];
-					struct PsiComplex u_jl_dag = psi_conjugate_complex(gate.data[tgt_j * gate_dim + l]);
+					struct PsiComplex u_jl_dag =
+							psi_conjugate_complex(gate.data[tgt_j * gate_dim + l]);
 					struct PsiComplex rho_kl = dm->data[src_i * dim + src_j];
 
-					sum = psi_add_complex(sum, psi_mul_complex(psi_mul_complex(u_ik, rho_kl), u_jl_dag));
+					sum = psi_add_complex(sum,
+					                      psi_mul_complex(psi_mul_complex(u_ik, rho_kl), u_jl_dag));
 				}
 
 			new_data[i * dim + j] = sum;
@@ -328,14 +350,15 @@ void psi_apply_unitary_density_matrix(struct PsiDensityMatrix *dm, struct PsiMat
 	dm->data = new_data;
 }
 
-void psi_apply_noise_channel(struct PsiDensityMatrix *dm, struct PsiNoiseChannel channel, size_t target)
+void psi_apply_noise_channel(struct PsiDensityMatrix* dm, struct PsiNoiseChannel channel,
+                             size_t target)
 {
 	assert(channel.num_qubits == 1);
 
 	size_t dim = dm->dim;
 	size_t target_bit = dm->num_qubits - 1 - target;
 
-	struct PsiComplex *new_data = calloc(dim * dim, sizeof(struct PsiComplex));
+	struct PsiComplex* new_data = calloc(dim * dim, sizeof(struct PsiComplex));
 	assert(new_data != NULL);
 
 	for (size_t op = 0; op < channel.operator_count; op++)
@@ -355,10 +378,12 @@ void psi_apply_noise_channel(struct PsiDensityMatrix *dm, struct PsiNoiseChannel
 						size_t src_j = (j & ~((size_t)1 << target_bit)) | (kj << target_bit);
 
 						struct PsiComplex k_elem = k.data[i_target * 2 + ki];
-						struct PsiComplex k_dag_elem = psi_conjugate_complex(k.data[j_target * 2 + kj]);
+						struct PsiComplex k_dag_elem =
+								psi_conjugate_complex(k.data[j_target * 2 + kj]);
 						struct PsiComplex rho_elem = dm->data[src_i * dim + src_j];
 
-						struct PsiComplex term = psi_mul_complex(psi_mul_complex(k_elem, rho_elem), k_dag_elem);
+						struct PsiComplex term =
+								psi_mul_complex(psi_mul_complex(k_elem, rho_elem), k_dag_elem);
 						new_data[i * dim + j] = psi_add_complex(new_data[i * dim + j], term);
 					}
 			}
@@ -368,7 +393,8 @@ void psi_apply_noise_channel(struct PsiDensityMatrix *dm, struct PsiNoiseChannel
 	dm->data = new_data;
 }
 
-double psi_measure_probability_density_matrix(struct PsiDensityMatrix dm, size_t qubit, size_t outcome)
+double psi_measure_probability_density_matrix(struct PsiDensityMatrix dm, size_t qubit,
+                                              size_t outcome)
 {
 	size_t target_bit = dm.num_qubits - 1 - qubit;
 	double prob = 0.0;
@@ -380,13 +406,16 @@ double psi_measure_probability_density_matrix(struct PsiDensityMatrix dm, size_t
 	return prob;
 }
 
-double psi_fidelity_density_matrix(struct PsiDensityMatrix dm, const struct PsiComplex *state)
+double psi_fidelity_density_matrix(struct PsiDensityMatrix dm, const struct PsiComplex* state)
 {
 	struct PsiComplex sum = psi_new_complex(0.0, 0.0);
 
 	for (size_t i = 0; i < dm.dim; i++)
 		for (size_t j = 0; j < dm.dim; j++)
-			sum = psi_add_complex(sum, psi_mul_complex(psi_mul_complex(psi_conjugate_complex(state[i]), dm.data[i * dm.dim + j]), state[j]));
+			sum = psi_add_complex(sum,
+			                      psi_mul_complex(psi_mul_complex(psi_conjugate_complex(state[i]),
+			                                                      dm.data[i * dm.dim + j]),
+			                                      state[j]));
 
 	return sum.real;
 }
diff --git a/src/core/quantum_components.c b/src/core/quantum_components.c
index 31d60dc..a62a911 100644
--- a/src/core/quantum_components.c
+++ b/src/core/quantum_components.c
@@ -14,21 +14,21 @@ struct PsiVector psi_new_state_1(void)
 	return psi_column_vector(psi_new_complex(0.0, 0.0), psi_new_complex(1.0, 0.0));
 }
 
-struct PsiQuantumGate psi_new_quantum_gate(const char *name, struct PsiMatrix matrix, size_t num_qubits)
+struct PsiQuantumGate psi_new_quantum_gate(const char* name, struct PsiMatrix matrix,
+                                           size_t num_qubits)
 {
 	size_t expected_dim = (size_t)1 << num_qubits;
 	assert(matrix.rows == expected_dim);
 	assert(matrix.cols == expected_dim);
 
-	return (struct PsiQuantumGate)
-	{
+	return (struct PsiQuantumGate){
 		name,
 		matrix,
 		num_qubits,
 	};
 }
 
-struct PsiQuantumGate psi_new_quantum_gate_from_matrix(const char *name, struct PsiMatrix matrix)
+struct PsiQuantumGate psi_new_quantum_gate_from_matrix(const char* name, struct PsiMatrix matrix)
 {
 	assert(matrix.rows == matrix.cols);
 
@@ -39,34 +39,32 @@ struct PsiQuantumGate psi_new_quantum_gate_from_matrix(const char *name, struct
 	while (((size_t)1 << num_qubits) < dim)
 		num_qubits++;
 
-	return (struct PsiQuantumGate)
-	{
+	return (struct PsiQuantumGate){
 		name,
 		matrix,
 		num_qubits,
 	};
 }
 
-void psi_free_quantum_gate(struct PsiQuantumGate *gate)
+void psi_free_quantum_gate(struct PsiQuantumGate* gate)
 {
 	psi_free_matrix(&gate->matrix);
 }
 
-struct PsiQuantumBit psi_new_quantum_bit(const char *name, struct PsiVector state)
+struct PsiQuantumBit psi_new_quantum_bit(const char* name, struct PsiVector state)
 {
-	return (struct PsiQuantumBit)
-	{
+	return (struct PsiQuantumBit){
 		name,
 		state,
 	};
 }
 
-void psi_free_quantum_bit(struct PsiQuantumBit *bit)
+void psi_free_quantum_bit(struct PsiQuantumBit* bit)
 {
 	psi_free_vector(&bit->state);
 }
 
-static void update_register(struct PsiQuantumRegister *reg)
+static void update_register(struct PsiQuantumRegister* reg)
 {
 	psi_free_vector(&reg->state_vector);
 
@@ -90,16 +88,16 @@ static void update_register(struct PsiQuantumRegister *reg)
 	psi_free_matrix(&result);
 }
 
-struct PsiQuantumRegister psi_new_quantum_register(const char *name, const char **names, size_t count)
+struct PsiQuantumRegister psi_new_quantum_register(const char* name, const char** names,
+                                                   size_t count)
 {
-	struct PsiQuantumBit *qubits = malloc(count * sizeof(struct PsiQuantumBit));
+	struct PsiQuantumBit* qubits = malloc(count * sizeof(struct PsiQuantumBit));
 	assert(qubits != NULL || count == 0);
 
 	for (size_t i = 0; i < count; i++)
 		qubits[i] = psi_new_quantum_bit(names[i], psi_new_state_0());
 
-	struct PsiQuantumRegister reg =
-	{
+	struct PsiQuantumRegister reg = {
 		name,
 		psi_new_vector(0, PSI_COLUMN_VECTOR),
 		qubits,
@@ -110,16 +108,16 @@ struct PsiQuantumRegister psi_new_quantum_register(const char *name, const char
 	return reg;
 }
 
-struct PsiQuantumRegister psi_new_quantum_register_from(const char *name, const struct PsiQuantumBit *bits, size_t count)
+struct PsiQuantumRegister
+psi_new_quantum_register_from(const char* name, const struct PsiQuantumBit* bits, size_t count)
 {
-	struct PsiQuantumBit *qubits = malloc(count * sizeof(struct PsiQuantumBit));
+	struct PsiQuantumBit* qubits = malloc(count * sizeof(struct PsiQuantumBit));
 	assert(qubits != NULL || count == 0);
 
 	for (size_t i = 0; i < count; i++)
 		qubits[i] = psi_new_quantum_bit(bits[i].name, psi_clone_vector(bits[i].state));
 
-	struct PsiQuantumRegister reg =
-	{
+	struct PsiQuantumRegister reg = {
 		name,
 		psi_new_vector(0, PSI_COLUMN_VECTOR),
 		qubits,
@@ -130,7 +128,7 @@ struct PsiQuantumRegister psi_new_quantum_register_from(const char *name, const
 	return reg;
 }
 
-void psi_free_quantum_register(struct PsiQuantumRegister *reg)
+void psi_free_quantum_register(struct PsiQuantumRegister* reg)
 {
 	for (size_t i = 0; i < reg->num_qubits; i++)
 		psi_free_quantum_bit(&reg->qubits[i]);
@@ -141,7 +139,7 @@ void psi_free_quantum_register(struct PsiQuantumRegister *reg)
 	psi_free_vector(&reg->state_vector);
 }
 
-static bool targets_contain(const size_t *targets, size_t count, size_t value)
+static bool targets_contain(const size_t* targets, size_t count, size_t value)
 {
 	for (size_t i = 0; i < count; i++)
 		if (targets[i] == value)
@@ -150,7 +148,8 @@ static bool targets_contain(const size_t *targets, size_t count, size_t value)
 	return false;
 }
 
-static struct PsiMatrix build_contiguous_operator(struct PsiQuantumRegister reg, struct PsiQuantumGate gate, size_t start_idx)
+static struct PsiMatrix build_contiguous_operator(struct PsiQuantumRegister reg,
+                                                  struct PsiQuantumGate gate, size_t start_idx)
 {
 	size_t n = reg.num_qubits;
 	size_t g = gate.num_qubits;
@@ -163,9 +162,8 @@ static struct PsiMatrix build_contiguous_operator(struct PsiQuantumRegister reg,
 		if (i > start_idx && i < start_idx + g)
 			continue;
 
-		struct PsiMatrix part = i == start_idx
-			? psi_clone_matrix(gate.matrix)
-			: psi_identity_matrix(2);
+		struct PsiMatrix part =
+				i == start_idx ? psi_clone_matrix(gate.matrix) : psi_identity_matrix(2);
 
 		if (!has_result)
 		{
@@ -186,7 +184,9 @@ static struct PsiMatrix build_contiguous_operator(struct PsiQuantumRegister reg,
 	return result;
 }
 
-static struct PsiMatrix build_full_operator(struct PsiQuantumRegister reg, struct PsiQuantumGate gate, const size_t *targets, size_t target_count)
+static struct PsiMatrix build_full_operator(struct PsiQuantumRegister reg,
+                                            struct PsiQuantumGate gate, const size_t* targets,
+                                            size_t target_count)
 {
 	size_t n = reg.num_qubits;
 	size_t g = gate.num_qubits;
@@ -238,13 +238,15 @@ static struct PsiMatrix build_full_operator(struct PsiQuantumRegister reg, struc
 			}
 
 			if (non_target_match)
-				result.data[row * result.cols + col] = psi_get_matrix(gate.matrix, target_row_bits, target_col_bits);
+				result.data[row * result.cols + col] =
+						psi_get_matrix(gate.matrix, target_row_bits, target_col_bits);
 		}
 
 	return result;
 }
 
-void psi_apply_gate(struct PsiQuantumRegister *reg, struct PsiQuantumGate gate, const size_t *targets, size_t target_count)
+void psi_apply_gate(struct PsiQuantumRegister* reg, struct PsiQuantumGate gate,
+                    const size_t* targets, size_t target_count)
 {
 	size_t n = reg->num_qubits;
 
diff --git a/src/core/runtime.c b/src/core/runtime.c
index 950ec7c..24cbe26 100644
--- a/src/core/runtime.c
+++ b/src/core/runtime.c
@@ -34,14 +34,9 @@ struct PsiRuntimeConfig psi_runtime_to_config(enum PsiRuntime runtime)
 
 	switch (runtime)
 	{
-		case PSI_RUNTIME_BASIC:
-			break;
-		case PSI_RUNTIME_BASIC_MT:
-			config.parallel = true;
-			break;
-		case PSI_RUNTIME_BATCHED:
-			config.batched = true;
-			break;
+		case PSI_RUNTIME_BASIC: break;
+		case PSI_RUNTIME_BASIC_MT: config.parallel = true; break;
+		case PSI_RUNTIME_BATCHED: config.batched = true; break;
 		case PSI_RUNTIME_BATCHED_MT:
 			config.batched = true;
 			config.parallel = true;
@@ -69,45 +64,126 @@ struct PsiRuntimeConfig psi_runtime_to_config(enum PsiRuntime runtime)
 	return config;
 }
 
-static bool op_to_kernel(struct PsiGateOp op, struct PsiKernel *out)
+static bool op_to_kernel(struct PsiGateOp op, struct PsiKernel* out)
 {
 	struct PsiMatrix matrix;
-	const char *name;
+	const char* name;
 
 	switch (op.kind)
 	{
-		case PSI_GATE_H: matrix = psi_hadamard_gate().matrix; name = "H"; break;
-		case PSI_GATE_X: matrix = psi_pauli_x_gate().matrix; name = "X"; break;
-		case PSI_GATE_Y: matrix = psi_pauli_y_gate().matrix; name = "Y"; break;
-		case PSI_GATE_Z: matrix = psi_pauli_z_gate().matrix; name = "Z"; break;
-		case PSI_GATE_S: matrix = psi_s_gate().matrix; name = "S"; break;
-		case PSI_GATE_T: matrix = psi_t_gate().matrix; name = "T"; break;
-		case PSI_GATE_SDG: matrix = psi_sdg_gate().matrix; name = "Sdg"; break;
-		case PSI_GATE_TDG: matrix = psi_tdg_gate().matrix; name = "Tdg"; break;
-		case PSI_GATE_SX: matrix = psi_sx_gate().matrix; name = "Sx"; break;
-		case PSI_GATE_SXDG: matrix = psi_sxdg_gate().matrix; name = "Sxdg"; break;
-		case PSI_GATE_RX: matrix = psi_rx_matrix(op.params[0]); name = "Rx"; break;
-		case PSI_GATE_RY: matrix = psi_ry_matrix(op.params[0]); name = "Ry"; break;
-		case PSI_GATE_RZ: matrix = psi_rz_matrix(op.params[0]); name = "Rz"; break;
-		case PSI_GATE_P: matrix = psi_p_matrix(op.params[0]); name = "P"; break;
-		case PSI_GATE_U1: matrix = psi_u1_matrix(op.params[0]); name = "U1"; break;
-		case PSI_GATE_U2: matrix = psi_u2_matrix(op.params[0], op.params[1]); name = "U2"; break;
-		case PSI_GATE_U3: matrix = psi_u3_matrix(op.params[0], op.params[1], op.params[2]); name = "U3"; break;
-		case PSI_GATE_CNOT: matrix = psi_cnot_gate().matrix; name = "CNOT"; break;
-		case PSI_GATE_CZ: matrix = psi_cz_gate().matrix; name = "CZ"; break;
-		case PSI_GATE_SWAP: matrix = psi_swap_gate().matrix; name = "SWAP"; break;
-		case PSI_GATE_CRX: matrix = psi_crx_matrix(op.params[0]); name = "CRx"; break;
-		case PSI_GATE_CRY: matrix = psi_cry_matrix(op.params[0]); name = "CRy"; break;
-		case PSI_GATE_CRZ: matrix = psi_crz_matrix(op.params[0]); name = "CRz"; break;
-		case PSI_GATE_CP: matrix = psi_cp_matrix(op.params[0]); name = "CP"; break;
-		case PSI_GATE_CCNOT: matrix = psi_toffoli_gate().matrix; name = "CCNOT"; break;
-		case PSI_GATE_CSWAP: matrix = psi_fredkin_gate().matrix; name = "CSWAP"; break;
+		case PSI_GATE_H:
+			matrix = psi_hadamard_gate().matrix;
+			name = "H";
+			break;
+		case PSI_GATE_X:
+			matrix = psi_pauli_x_gate().matrix;
+			name = "X";
+			break;
+		case PSI_GATE_Y:
+			matrix = psi_pauli_y_gate().matrix;
+			name = "Y";
+			break;
+		case PSI_GATE_Z:
+			matrix = psi_pauli_z_gate().matrix;
+			name = "Z";
+			break;
+		case PSI_GATE_S:
+			matrix = psi_s_gate().matrix;
+			name = "S";
+			break;
+		case PSI_GATE_T:
+			matrix = psi_t_gate().matrix;
+			name = "T";
+			break;
+		case PSI_GATE_SDG:
+			matrix = psi_sdg_gate().matrix;
+			name = "Sdg";
+			break;
+		case PSI_GATE_TDG:
+			matrix = psi_tdg_gate().matrix;
+			name = "Tdg";
+			break;
+		case PSI_GATE_SX:
+			matrix = psi_sx_gate().matrix;
+			name = "Sx";
+			break;
+		case PSI_GATE_SXDG:
+			matrix = psi_sxdg_gate().matrix;
+			name = "Sxdg";
+			break;
+		case PSI_GATE_RX:
+			matrix = psi_rx_matrix(op.params[0]);
+			name = "Rx";
+			break;
+		case PSI_GATE_RY:
+			matrix = psi_ry_matrix(op.params[0]);
+			name = "Ry";
+			break;
+		case PSI_GATE_RZ:
+			matrix = psi_rz_matrix(op.params[0]);
+			name = "Rz";
+			break;
+		case PSI_GATE_P:
+			matrix = psi_p_matrix(op.params[0]);
+			name = "P";
+			break;
+		case PSI_GATE_U1:
+			matrix = psi_u1_matrix(op.params[0]);
+			name = "U1";
+			break;
+		case PSI_GATE_U2:
+			matrix = psi_u2_matrix(op.params[0], op.params[1]);
+			name = "U2";
+			break;
+		case PSI_GATE_U3:
+			matrix = psi_u3_matrix(op.params[0], op.params[1], op.params[2]);
+			name = "U3";
+			break;
+		case PSI_GATE_CNOT:
+			matrix = psi_cnot_gate().matrix;
+			name = "CNOT";
+			break;
+		case PSI_GATE_CZ:
+			matrix = psi_cz_gate().matrix;
+			name = "CZ";
+			break;
+		case PSI_GATE_SWAP:
+			matrix = psi_swap_gate().matrix;
+			name = "SWAP";
+			break;
+		case PSI_GATE_CRX:
+			matrix = psi_crx_matrix(op.params[0]);
+			name = "CRx";
+			break;
+		case PSI_GATE_CRY:
+			matrix = psi_cry_matrix(op.params[0]);
+			name = "CRy";
+			break;
+		case PSI_GATE_CRZ:
+			matrix = psi_crz_matrix(op.params[0]);
+			name = "CRz";
+			break;
+		case PSI_GATE_CP:
+			matrix = psi_cp_matrix(op.params[0]);
+			name = "CP";
+			break;
+		case PSI_GATE_CCNOT:
+			matrix = psi_toffoli_gate().matrix;
+			name = "CCNOT";
+			break;
+		case PSI_GATE_CSWAP:
+			matrix = psi_fredkin_gate().matrix;
+			name = "CSWAP";
+			break;
 		case PSI_GATE_MEASURE: return false;
-		case PSI_GATE_CUSTOM: matrix = psi_to_quantum_gate(*op.custom).matrix; name = "Custom"; break;
+		case PSI_GATE_CUSTOM:
+			matrix = psi_to_quantum_gate(*op.custom).matrix;
+			name = "Custom";
+			break;
 	}
 
 	size_t target_count;
-	const size_t *targets = psi_gate_op_quantum_targets(&op, &target_count);
+	const size_t* targets = psi_gate_op_quantum_targets(&op, &target_count);
 	*out = psi_new_kernel(name, matrix, targets, target_count);
 
 	return true;
@@ -122,7 +198,8 @@ static struct PsiVector new_zero_state(size_t num_qubits)
 	return state;
 }
 
-static void execute_kernels(struct PsiVector *state, const struct PsiKernel *kernels, size_t count, size_t num_qubits, struct PsiRuntimeConfig config)
+static void execute_kernels(struct PsiVector* state, const struct PsiKernel* kernels, size_t count,
+                            size_t num_qubits, struct PsiRuntimeConfig config)
 {
 	bool use_parallel = config.parallel && num_qubits >= config.parallel_threshold;
 
@@ -138,7 +215,8 @@ static void execute_kernels(struct PsiVector *state, const struct PsiKernel *ker
 			};
 
 			if (use_parallel)
-				psi_apply_single_qubit_gate_simd_parallel(state->data, gate, kernel.targets[0], num_qubits);
+				psi_apply_single_qubit_gate_simd_parallel(state->data, gate, kernel.targets[0],
+				                                          num_qubits);
 			else
 				psi_apply_single_qubit_gate_simd(state->data, gate, kernel.targets[0], num_qubits);
 		}
@@ -147,7 +225,8 @@ static void execute_kernels(struct PsiVector *state, const struct PsiKernel *ker
 	}
 }
 
-struct PsiVector psi_compute_runtime_config(struct PsiRuntimeConfig config, size_t num_qubits, const struct PsiGateOp *operations, size_t op_count)
+struct PsiVector psi_compute_runtime_config(struct PsiRuntimeConfig config, size_t num_qubits,
+                                            const struct PsiGateOp* operations, size_t op_count)
 {
 	struct PsiVector state = new_zero_state(num_qubits);
 
@@ -185,42 +264,47 @@ struct PsiVector psi_compute_runtime_config(struct PsiRuntimeConfig config, size
 	return state;
 }
 
-struct PsiVector psi_compute_runtime(enum PsiRuntime runtime, size_t num_qubits, const struct PsiGateOp *operations, size_t op_count)
+struct PsiVector psi_compute_runtime(enum PsiRuntime runtime, size_t num_qubits,
+                                     const struct PsiGateOp* operations, size_t op_count)
 {
-	return psi_compute_runtime_config(psi_runtime_to_config(runtime), num_qubits, operations, op_count);
+	return psi_compute_runtime_config(psi_runtime_to_config(runtime), num_qubits, operations,
+	                                  op_count);
 }
 
-const struct PsiVector *psi_compute_circuit_with_config(struct PsiQuantumCircuit *circuit, struct PsiRuntimeConfig config)
+const struct PsiVector* psi_compute_circuit_with_config(struct PsiQuantumCircuit* circuit,
+                                                        struct PsiRuntimeConfig config)
 {
 	if (!circuit->is_computed)
 	{
 		psi_free_vector(&circuit->computed_state);
-		circuit->computed_state = psi_compute_runtime_config(config, circuit->num_qubits, circuit->operations, circuit->operation_count);
+		circuit->computed_state = psi_compute_runtime_config(
+				config, circuit->num_qubits, circuit->operations, circuit->operation_count);
 		circuit->is_computed = true;
 	}
 
 	return &circuit->computed_state;
 }
 
-const struct PsiVector *psi_compute_circuit_with(struct PsiQuantumCircuit *circuit, enum PsiRuntime runtime)
+const struct PsiVector* psi_compute_circuit_with(struct PsiQuantumCircuit* circuit,
+                                                 enum PsiRuntime runtime)
 {
 	return psi_compute_circuit_with_config(circuit, psi_runtime_to_config(runtime));
 }
 
-const struct PsiVector *psi_compute_circuit(struct PsiQuantumCircuit *circuit)
+const struct PsiVector* psi_compute_circuit(struct PsiQuantumCircuit* circuit)
 {
 	return psi_compute_circuit_with(circuit, PSI_RUNTIME_BASIC);
 }
 
-double psi_circuit_probability(struct PsiQuantumCircuit *circuit, size_t state_index)
+double psi_circuit_probability(struct PsiQuantumCircuit* circuit, size_t state_index)
 {
-	const struct PsiVector *state = psi_compute_circuit(circuit);
+	const struct PsiVector* state = psi_compute_circuit(circuit);
 	return psi_norm2_complex(state->data[state_index]);
 }
 
-void psi_circuit_probabilities(struct PsiQuantumCircuit *circuit, double *out)
+void psi_circuit_probabilities(struct PsiQuantumCircuit* circuit, double* out)
 {
-	const struct PsiVector *state = psi_compute_circuit(circuit);
+	const struct PsiVector* state = psi_compute_circuit(circuit);
 	size_t dim = (size_t)1 << circuit->num_qubits;
 	for (size_t i = 0; i < dim; i++)
 		out[i] = psi_norm2_complex(state->data[i]);
diff --git a/src/maths/complex.c b/src/maths/complex.c
index 3115fdd..44e255f 100644
--- a/src/maths/complex.c
+++ b/src/maths/complex.c
@@ -4,8 +4,7 @@
 
 struct PsiComplex psi_new_complex(double real, double imaginary)
 {
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		real,
 		imaginary,
 	};
@@ -13,8 +12,7 @@ struct PsiComplex psi_new_complex(double real, double imaginary)
 
 struct PsiComplex psi_new_complex_from_real(double real)
 {
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		real,
 		0.0,
 	};
@@ -22,8 +20,7 @@ struct PsiComplex psi_new_complex_from_real(double real)
 
 struct PsiComplex psi_conjugate_complex(struct PsiComplex z)
 {
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		z.real,
 		-z.imaginary,
 	};
@@ -31,8 +28,7 @@ struct PsiComplex psi_conjugate_complex(struct PsiComplex z)
 
 struct PsiComplex psi_neg_complex(struct PsiComplex z)
 {
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		-z.real,
 		-z.imaginary,
 	};
@@ -59,8 +55,7 @@ struct PsiComplex psi_sqrt_complex(struct PsiComplex z)
 	double half_theta = psi_phase_complex(z) / 2.0;
 	double sqrt_r = sqrt(r);
 
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		sqrt_r * cos(half_theta),
 		sqrt_r * sin(half_theta),
 	};
@@ -68,8 +63,7 @@ struct PsiComplex psi_sqrt_complex(struct PsiComplex z)
 
 struct PsiComplex psi_add_complex(struct PsiComplex a, struct PsiComplex b)
 {
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		a.real + b.real,
 		a.imaginary + b.imaginary,
 	};
@@ -77,8 +71,7 @@ struct PsiComplex psi_add_complex(struct PsiComplex a, struct PsiComplex b)
 
 struct PsiComplex psi_sub_complex(struct PsiComplex a, struct PsiComplex b)
 {
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		a.real - b.real,
 		a.imaginary - b.imaginary,
 	};
@@ -87,8 +80,7 @@ struct PsiComplex psi_sub_complex(struct PsiComplex a, struct PsiComplex b)
 struct PsiComplex psi_mul_complex(struct PsiComplex a, struct PsiComplex b)
 {
 	// (a + bi)(c + di) = (ac - bd) + (ad + bc)i
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		a.real * b.real - a.imaginary * b.imaginary,
 		a.real * b.imaginary + a.imaginary * b.real,
 	};
@@ -98,8 +90,7 @@ struct PsiComplex psi_div_complex(struct PsiComplex a, struct PsiComplex b)
 {
 	// (a + bi) / (c + di) = ((ac + bd) + (bc - ad)i) / (c² + d²)
 	double denom = b.real * b.real + b.imaginary * b.imaginary;
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		(a.real * b.real + a.imaginary * b.imaginary) / denom,
 		(a.imaginary * b.real - a.real * b.imaginary) / denom,
 	};
@@ -107,8 +98,7 @@ struct PsiComplex psi_div_complex(struct PsiComplex a, struct PsiComplex b)
 
 struct PsiComplex psi_add_complex_real(struct PsiComplex a, double b)
 {
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		a.real + b,
 		a.imaginary,
 	};
@@ -116,8 +106,7 @@ struct PsiComplex psi_add_complex_real(struct PsiComplex a, double b)
 
 struct PsiComplex psi_sub_complex_real(struct PsiComplex a, double b)
 {
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		a.real - b,
 		a.imaginary,
 	};
@@ -125,8 +114,7 @@ struct PsiComplex psi_sub_complex_real(struct PsiComplex a, double b)
 
 struct PsiComplex psi_mul_complex_real(struct PsiComplex a, double b)
 {
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		a.real * b,
 		a.imaginary * b,
 	};
@@ -134,8 +122,7 @@ struct PsiComplex psi_mul_complex_real(struct PsiComplex a, double b)
 
 struct PsiComplex psi_div_complex_real(struct PsiComplex a, double b)
 {
-	return (struct PsiComplex)
-	{
+	return (struct PsiComplex){
 		a.real / b,
 		a.imaginary / b,
 	};
diff --git a/src/maths/format.c b/src/maths/format.c
index ad80ff0..f3db83d 100644
--- a/src/maths/format.c
+++ b/src/maths/format.c
@@ -15,10 +15,10 @@ static bool approx_eq(double a, double b)
 	return fabs(a - b) < EPSILON;
 }
 
-static bool real_symbolic(double v, char *out, size_t cap)
+static bool real_symbolic(double v, char* out, size_t cap)
 {
 	double abs_v = fabs(v);
-	const char *sign = v < 0.0 ? "-" : "";
+	const char* sign = v < 0.0 ? "-" : "";
 
 	if (approx_eq(abs_v, 0.0))
 	{
@@ -26,7 +26,7 @@ static bool real_symbolic(double v, char *out, size_t cap)
 		return true;
 	}
 
-	const char *sym = NULL;
+	const char* sym = NULL;
 	if (approx_eq(abs_v, 1.0))
 		sym = "1";
 	else if (approx_eq(abs_v, 0.5))
@@ -58,7 +58,7 @@ static bool real_symbolic(double v, char *out, size_t cap)
 	return true;
 }
 
-char *psi_format_amplitude(struct PsiComplex c, char *out, size_t cap)
+char* psi_format_amplitude(struct PsiComplex c, char* out, size_t cap)
 {
 	double re = c.real;
 	double im = c.imaginary;
@@ -108,7 +108,7 @@ char *psi_format_amplitude(struct PsiComplex c, char *out, size_t cap)
 	}
 	else
 	{
-		const char *sign = im > 0.0 ? "+" : "-";
+		const char* sign = im > 0.0 ? "+" : "-";
 		char sym[32];
 		if (real_symbolic(fabs(im), sym, sizeof sym))
 			snprintf(im_str, sizeof im_str, "%s%si", sign, sym);
@@ -120,7 +120,7 @@ char *psi_format_amplitude(struct PsiComplex c, char *out, size_t cap)
 	return out;
 }
 
-char *psi_format_probability(double p, char *out, size_t cap)
+char* psi_format_probability(double p, char* out, size_t cap)
 {
 	if (approx_eq(p, 0.0))
 		snprintf(out, cap, "%s", "0");
diff --git a/src/maths/matrix.c b/src/maths/matrix.c
index 44b953c..2415395 100644
--- a/src/maths/matrix.c
+++ b/src/maths/matrix.c
@@ -6,18 +6,17 @@
 
 struct PsiMatrix psi_new_matrix(size_t rows, size_t cols)
 {
-	struct PsiComplex *data = calloc(rows * cols, sizeof(struct PsiComplex));
+	struct PsiComplex* data = calloc(rows * cols, sizeof(struct PsiComplex));
 	assert(data != NULL || rows * cols == 0);
 
-	return (struct PsiMatrix)
-	{
+	return (struct PsiMatrix){
 		data,
 		rows,
 		cols,
 	};
 }
 
-struct PsiMatrix psi_new_matrix_from(const struct PsiComplex *data, size_t rows, size_t cols)
+struct PsiMatrix psi_new_matrix_from(const struct PsiComplex* data, size_t rows, size_t cols)
 {
 	struct PsiMatrix m = psi_new_matrix(rows, cols);
 	memcpy(m.data, data, rows * cols * sizeof(struct PsiComplex));
@@ -39,7 +38,7 @@ struct PsiMatrix psi_identity_matrix(size_t size)
 	return m;
 }
 
-void psi_free_matrix(struct PsiMatrix *m)
+void psi_free_matrix(struct PsiMatrix* m)
 {
 	free(m->data);
 	m->data = NULL;
@@ -53,7 +52,7 @@ struct PsiComplex psi_get_matrix(struct PsiMatrix m, size_t row, size_t col)
 	return m.data[row * m.cols + col];
 }
 
-void psi_set_matrix(struct PsiMatrix *m, size_t row, size_t col, struct PsiComplex value)
+void psi_set_matrix(struct PsiMatrix* m, size_t row, size_t col, struct PsiComplex value)
 {
 	assert(row < m->rows && col < m->cols);
 	m->data[row * m->cols + col] = value;
@@ -69,7 +68,8 @@ struct PsiMatrix psi_dot_matrix(struct PsiMatrix a, struct PsiMatrix b)
 		{
 			struct PsiComplex sum = psi_new_complex(0.0, 0.0);
 			for (size_t k = 0; k < a.cols; k++)
-				sum = psi_add_complex(sum, psi_mul_complex(a.data[i * a.cols + k], b.data[k * b.cols + j]));
+				sum = psi_add_complex(
+						sum, psi_mul_complex(a.data[i * a.cols + k], b.data[k * b.cols + j]));
 
 			result.data[i * result.cols + j] = sum;
 		}
@@ -89,7 +89,8 @@ struct PsiMatrix psi_kronecker_matrix(struct PsiMatrix a, struct PsiMatrix b)
 				{
 					size_t row = i * b.rows + k;
 					size_t col = j * b.cols + l;
-					result.data[row * result.cols + col] = psi_mul_complex(a_val, b.data[k * b.cols + l]);
+					result.data[row * result.cols + col] =
+							psi_mul_complex(a_val, b.data[k * b.cols + l]);
 				}
 		}
 
diff --git a/src/maths/simd.c b/src/maths/simd.c
new file mode 100644
index 0000000..9e83dcc
--- /dev/null
+++ b/src/maths/simd.c
@@ -0,0 +1,324 @@
+#include "maths/simd.h"
+
+#include <assert.h>
+#include <stdlib.h>
+
+#if defined(__x86_64__) || defined(__i386__)
+#include <immintrin.h>
+#elif defined(__aarch64__)
+#include <arm_neon.h>
+#endif
+
+enum PsiSimdCapability psi_detect_simd(void)
+{
+#if defined(__x86_64__) || defined(__i386__)
+	if (__builtin_cpu_supports("avx512f") && __builtin_cpu_supports("avx512dq"))
+		return PSI_SIMD_AVX512;
+	if (__builtin_cpu_supports("avx2") && __builtin_cpu_supports("fma"))
+		return PSI_SIMD_AVX2;
+
+	return PSI_SIMD_NONE;
+#elif defined(__aarch64__)
+	return PSI_SIMD_NEON;
+#else
+	return PSI_SIMD_NONE;
+#endif
+}
+
+const char* psi_simd_name(enum PsiSimdCapability cap)
+{
+	switch (cap)
+	{
+		case PSI_SIMD_NONE: return "Scalar";
+		case PSI_SIMD_AVX2: return "AVX2+FMA";
+		case PSI_SIMD_AVX512: return "AVX-512";
+		case PSI_SIMD_NEON: return "NEON";
+	}
+
+	return "Scalar";
+}
+
+static void apply_pair(struct PsiComplex* state, size_t i, size_t j, struct PsiComplex g00,
+                       struct PsiComplex g01, struct PsiComplex g10, struct PsiComplex g11)
+{
+	struct PsiComplex s0 = state[i];
+	struct PsiComplex s1 = state[j];
+
+	state[i] = psi_new_complex(s0.real * g00.real - s0.imaginary * g00.imaginary +
+	                                   s1.real * g01.real - s1.imaginary * g01.imaginary,
+	                           s0.real * g00.imaginary + s0.imaginary * g00.real +
+	                                   s1.real * g01.imaginary + s1.imaginary * g01.real);
+	state[j] = psi_new_complex(s0.real * g10.real - s0.imaginary * g10.imaginary +
+	                                   s1.real * g11.real - s1.imaginary * g11.imaginary,
+	                           s0.real * g10.imaginary + s0.imaginary * g10.real +
+	                                   s1.real * g11.imaginary + s1.imaginary * g11.real);
+}
+
+static void apply_scalar(struct PsiComplex* state, const struct PsiComplex gate[2][2],
+                         size_t target, size_t num_qubits)
+{
+	size_t target_bit = num_qubits - 1 - target;
+	size_t step = (size_t)1 << target_bit;
+	size_t dim = (size_t)1 << num_qubits;
+
+	for (size_t i = 0; i < dim; i++)
+	{
+		if ((i >> target_bit) & 1)
+			continue;
+
+		apply_pair(state, i, i | step, gate[0][0], gate[0][1], gate[1][0], gate[1][1]);
+	}
+}
+
+static size_t (*build_pairs(size_t dim, size_t target_bit, size_t step, size_t* out_count))[2]
+{
+	size_t (*pairs)[2] = malloc((dim / 2) * sizeof(*pairs));
+	assert(pairs != NULL || dim == 0);
+
+	size_t n = 0;
+	for (size_t i = 0; i < dim; i++)
+		if (((i >> target_bit) & 1) == 0)
+		{
+			pairs[n][0] = i;
+			pairs[n][1] = i | step;
+			n++;
+		}
+
+	*out_count = n;
+	return pairs;
+}
+
+#if defined(__aarch64__)
+static void apply_neon(struct PsiComplex* state, const struct PsiComplex gate[2][2], size_t target,
+                       size_t num_qubits)
+{
+	size_t target_bit = num_qubits - 1 - target;
+	size_t step = (size_t)1 << target_bit;
+	size_t dim = (size_t)1 << num_qubits;
+
+	struct PsiComplex g00 = gate[0][0];
+	struct PsiComplex g01 = gate[0][1];
+	struct PsiComplex g10 = gate[1][0];
+	struct PsiComplex g11 = gate[1][1];
+
+	size_t np;
+	size_t (*pairs)[2] = build_pairs(dim, target_bit, step, &np);
+	size_t chunks = np / 2;
+
+	float64x2_t g00_re = vdupq_n_f64(g00.real);
+	float64x2_t g00_im = vdupq_n_f64(g00.imaginary);
+	float64x2_t g01_re = vdupq_n_f64(g01.real);
+	float64x2_t g01_im = vdupq_n_f64(g01.imaginary);
+	float64x2_t g10_re = vdupq_n_f64(g10.real);
+	float64x2_t g10_im = vdupq_n_f64(g10.imaginary);
+	float64x2_t g11_re = vdupq_n_f64(g11.real);
+	float64x2_t g11_im = vdupq_n_f64(g11.imaginary);
+
+	for (size_t chunk = 0; chunk < chunks; chunk++)
+	{
+		size_t i0 = pairs[chunk * 2][0];
+		size_t j0 = pairs[chunk * 2][1];
+		size_t i1 = pairs[chunk * 2 + 1][0];
+		size_t j1 = pairs[chunk * 2 + 1][1];
+
+		double s0r[2] = { state[i0].real, state[i1].real };
+		double s0i[2] = { state[i0].imaginary, state[i1].imaginary };
+		double s1r[2] = { state[j0].real, state[j1].real };
+		double s1i[2] = { state[j0].imaginary, state[j1].imaginary };
+
+		float64x2_t s0_re = vld1q_f64(s0r);
+		float64x2_t s0_im = vld1q_f64(s0i);
+		float64x2_t s1_re = vld1q_f64(s1r);
+		float64x2_t s1_im = vld1q_f64(s1i);
+
+		float64x2_t new0_re = vaddq_f64(vfmsq_f64(vmulq_f64(s0_re, g00_re), s0_im, g00_im),
+		                                vfmsq_f64(vmulq_f64(s1_re, g01_re), s1_im, g01_im));
+		float64x2_t new0_im = vaddq_f64(vfmaq_f64(vmulq_f64(s0_re, g00_im), s0_im, g00_re),
+		                                vfmaq_f64(vmulq_f64(s1_re, g01_im), s1_im, g01_re));
+		float64x2_t new1_re = vaddq_f64(vfmsq_f64(vmulq_f64(s0_re, g10_re), s0_im, g10_im),
+		                                vfmsq_f64(vmulq_f64(s1_re, g11_re), s1_im, g11_im));
+		float64x2_t new1_im = vaddq_f64(vfmaq_f64(vmulq_f64(s0_re, g10_im), s0_im, g10_re),
+		                                vfmaq_f64(vmulq_f64(s1_re, g11_im), s1_im, g11_re));
+
+		state[i0] = psi_new_complex(vgetq_lane_f64(new0_re, 0), vgetq_lane_f64(new0_im, 0));
+		state[j0] = psi_new_complex(vgetq_lane_f64(new1_re, 0), vgetq_lane_f64(new1_im, 0));
+		state[i1] = psi_new_complex(vgetq_lane_f64(new0_re, 1), vgetq_lane_f64(new0_im, 1));
+		state[j1] = psi_new_complex(vgetq_lane_f64(new1_re, 1), vgetq_lane_f64(new1_im, 1));
+	}
+
+	for (size_t p = chunks * 2; p < np; p++)
+		apply_pair(state, pairs[p][0], pairs[p][1], g00, g01, g10, g11);
+
+	free(pairs);
+}
+#endif
+
+#if defined(__x86_64__) || defined(__i386__)
+__attribute__((target("avx2,fma"))) static void apply_avx2(struct PsiComplex* state,
+                                                           const struct PsiComplex gate[2][2],
+                                                           size_t target, size_t num_qubits)
+{
+	size_t target_bit = num_qubits - 1 - target;
+	size_t step = (size_t)1 << target_bit;
+	size_t dim = (size_t)1 << num_qubits;
+
+	struct PsiComplex g00 = gate[0][0];
+	struct PsiComplex g01 = gate[0][1];
+	struct PsiComplex g10 = gate[1][0];
+	struct PsiComplex g11 = gate[1][1];
+
+	size_t np;
+	size_t (*pairs)[2] = build_pairs(dim, target_bit, step, &np);
+	size_t chunks = np / 2;
+
+	for (size_t chunk = 0; chunk < chunks; chunk++)
+	{
+		size_t i0 = pairs[chunk * 2][0];
+		size_t j0 = pairs[chunk * 2][1];
+		size_t i1 = pairs[chunk * 2 + 1][0];
+		size_t j1 = pairs[chunk * 2 + 1][1];
+
+		__m256d s0_re =
+				_mm256_set_pd(state[j1].real, state[i1].real, state[j0].real, state[i0].real);
+		__m256d s0_im = _mm256_set_pd(state[j1].imaginary, state[i1].imaginary, state[j0].imaginary,
+		                              state[i0].imaginary);
+
+		__m256d g_re_0 = _mm256_set_pd(g01.real, g00.real, g01.real, g00.real);
+		__m256d g_im_0 = _mm256_set_pd(g01.imaginary, g00.imaginary, g01.imaginary, g00.imaginary);
+		__m256d g_re_1 = _mm256_set_pd(g11.real, g10.real, g11.real, g10.real);
+		__m256d g_im_1 = _mm256_set_pd(g11.imaginary, g10.imaginary, g11.imaginary, g10.imaginary);
+
+		__m256d prod0_re = _mm256_fmsub_pd(s0_re, g_re_0, _mm256_mul_pd(s0_im, g_im_0));
+		__m256d prod0_im = _mm256_fmadd_pd(s0_re, g_im_0, _mm256_mul_pd(s0_im, g_re_0));
+		__m256d prod1_re = _mm256_fmsub_pd(s0_re, g_re_1, _mm256_mul_pd(s0_im, g_im_1));
+		__m256d prod1_im = _mm256_fmadd_pd(s0_re, g_im_1, _mm256_mul_pd(s0_im, g_re_1));
+
+		double res0_re[4];
+		double res0_im[4];
+		double res1_re[4];
+		double res1_im[4];
+
+		_mm256_storeu_pd(res0_re, prod0_re);
+		_mm256_storeu_pd(res0_im, prod0_im);
+		_mm256_storeu_pd(res1_re, prod1_re);
+		_mm256_storeu_pd(res1_im, prod1_im);
+
+		state[i0] = psi_new_complex(res0_re[0] + res0_re[1], res0_im[0] + res0_im[1]);
+		state[j0] = psi_new_complex(res1_re[0] + res1_re[1], res1_im[0] + res1_im[1]);
+		state[i1] = psi_new_complex(res0_re[2] + res0_re[3], res0_im[2] + res0_im[3]);
+		state[j1] = psi_new_complex(res1_re[2] + res1_re[3], res1_im[2] + res1_im[3]);
+	}
+
+	for (size_t p = chunks * 2; p < np; p++)
+		apply_pair(state, pairs[p][0], pairs[p][1], g00, g01, g10, g11);
+
+	free(pairs);
+}
+
+__attribute__((target("avx512f,avx512dq"))) static void
+apply_avx512(struct PsiComplex* state, const struct PsiComplex gate[2][2], size_t target,
+             size_t num_qubits)
+{
+	size_t target_bit = num_qubits - 1 - target;
+	size_t step = (size_t)1 << target_bit;
+	size_t dim = (size_t)1 << num_qubits;
+
+	struct PsiComplex g00 = gate[0][0];
+	struct PsiComplex g01 = gate[0][1];
+	struct PsiComplex g10 = gate[1][0];
+	struct PsiComplex g11 = gate[1][1];
+
+	size_t np;
+	size_t (*pairs)[2] = build_pairs(dim, target_bit, step, &np);
+	size_t chunks = np / 4;
+
+	for (size_t chunk = 0; chunk < chunks; chunk++)
+	{
+		size_t base = chunk * 4;
+		size_t i0 = pairs[base][0], j0 = pairs[base][1];
+		size_t i1 = pairs[base + 1][0], j1 = pairs[base + 1][1];
+		size_t i2 = pairs[base + 2][0], j2 = pairs[base + 2][1];
+		size_t i3 = pairs[base + 3][0], j3 = pairs[base + 3][1];
+
+		__m512d s0_re =
+				_mm512_set_pd(state[j3].real, state[i3].real, state[j2].real, state[i2].real,
+				              state[j1].real, state[i1].real, state[j0].real, state[i0].real);
+		__m512d s0_im = _mm512_set_pd(state[j3].imaginary, state[i3].imaginary, state[j2].imaginary,
+		                              state[i2].imaginary, state[j1].imaginary, state[i1].imaginary,
+		                              state[j0].imaginary, state[i0].imaginary);
+
+		__m512d g_re_0 = _mm512_set_pd(g01.real, g00.real, g01.real, g00.real, g01.real, g00.real,
+		                               g01.real, g00.real);
+		__m512d g_im_0 = _mm512_set_pd(g01.imaginary, g00.imaginary, g01.imaginary, g00.imaginary,
+		                               g01.imaginary, g00.imaginary, g01.imaginary, g00.imaginary);
+		__m512d g_re_1 = _mm512_set_pd(g11.real, g10.real, g11.real, g10.real, g11.real, g10.real,
+		                               g11.real, g10.real);
+		__m512d g_im_1 = _mm512_set_pd(g11.imaginary, g10.imaginary, g11.imaginary, g10.imaginary,
+		                               g11.imaginary, g10.imaginary, g11.imaginary, g10.imaginary);
+
+		__m512d prod0_re = _mm512_fmsub_pd(s0_re, g_re_0, _mm512_mul_pd(s0_im, g_im_0));
+		__m512d prod0_im = _mm512_fmadd_pd(s0_re, g_im_0, _mm512_mul_pd(s0_im, g_re_0));
+		__m512d prod1_re = _mm512_fmsub_pd(s0_re, g_re_1, _mm512_mul_pd(s0_im, g_im_1));
+		__m512d prod1_im = _mm512_fmadd_pd(s0_re, g_im_1, _mm512_mul_pd(s0_im, g_re_1));
+
+		double res0_re[8];
+		double res0_im[8];
+		double res1_re[8];
+		double res1_im[8];
+
+		_mm512_storeu_pd(res0_re, prod0_re);
+		_mm512_storeu_pd(res0_im, prod0_im);
+		_mm512_storeu_pd(res1_re, prod1_re);
+		_mm512_storeu_pd(res1_im, prod1_im);
+
+		state[i0] = psi_new_complex(res0_re[0] + res0_re[1], res0_im[0] + res0_im[1]);
+		state[j0] = psi_new_complex(res1_re[0] + res1_re[1], res1_im[0] + res1_im[1]);
+		state[i1] = psi_new_complex(res0_re[2] + res0_re[3], res0_im[2] + res0_im[3]);
+		state[j1] = psi_new_complex(res1_re[2] + res1_re[3], res1_im[2] + res1_im[3]);
+		state[i2] = psi_new_complex(res0_re[4] + res0_re[5], res0_im[4] + res0_im[5]);
+		state[j2] = psi_new_complex(res1_re[4] + res1_re[5], res1_im[4] + res1_im[5]);
+		state[i3] = psi_new_complex(res0_re[6] + res0_re[7], res0_im[6] + res0_im[7]);
+		state[j3] = psi_new_complex(res1_re[6] + res1_re[7], res1_im[6] + res1_im[7]);
+	}
+
+	for (size_t p = chunks * 4; p < np; p++)
+		apply_pair(state, pairs[p][0], pairs[p][1], g00, g01, g10, g11);
+
+	free(pairs);
+}
+#endif
+
+void psi_apply_single_qubit_gate_simd(struct PsiComplex* state, const struct PsiComplex gate[2][2],
+                                      size_t target, size_t num_qubits)
+{
+	enum PsiSimdCapability cap = psi_detect_simd();
+
+#if defined(__x86_64__) || defined(__i386__)
+	if (cap == PSI_SIMD_AVX512)
+	{
+		apply_avx512(state, gate, target, num_qubits);
+		return;
+	}
+	if (cap == PSI_SIMD_AVX2)
+	{
+		apply_avx2(state, gate, target, num_qubits);
+		return;
+	}
+#elif defined(__aarch64__)
+	if (cap == PSI_SIMD_NEON)
+	{
+		apply_neon(state, gate, target, num_qubits);
+		return;
+	}
+#endif
+
+	(void)cap;
+	apply_scalar(state, gate, target, num_qubits);
+}
+
+void psi_apply_single_qubit_gate_simd_parallel(struct PsiComplex* state,
+                                               const struct PsiComplex gate[2][2], size_t target,
+                                               size_t num_qubits)
+{
+	apply_scalar(state, gate, target, num_qubits);
+}
diff --git a/src/maths/vector.c b/src/maths/vector.c
index b9713ae..1c5c4f9 100644
--- a/src/maths/vector.c
+++ b/src/maths/vector.c
@@ -14,18 +14,18 @@ static enum PsiVectorKind flip_kind(enum PsiVectorKind kind)
 
 struct PsiVector psi_new_vector(size_t size, enum PsiVectorKind kind)
 {
-	struct PsiComplex *data = calloc(size, sizeof(struct PsiComplex));
+	struct PsiComplex* data = calloc(size, sizeof(struct PsiComplex));
 	assert(data != NULL || size == 0);
 
-	return (struct PsiVector)
-	{
+	return (struct PsiVector){
 		data,
 		size,
 		kind,
 	};
 }
 
-struct PsiVector psi_new_vector_from(const struct PsiComplex *data, size_t size, enum PsiVectorKind kind)
+struct PsiVector psi_new_vector_from(const struct PsiComplex* data, size_t size,
+                                     enum PsiVectorKind kind)
 {
 	struct PsiVector v = psi_new_vector(size, kind);
 	memcpy(v.data, data, size * sizeof(struct PsiComplex));
@@ -38,7 +38,7 @@ struct PsiVector psi_clone_vector(struct PsiVector v)
 	return psi_new_vector_from(v.data, v.size, v.kind);
 }
 
-void psi_free_vector(struct PsiVector *v)
+void psi_free_vector(struct PsiVector* v)
 {
 	free(v->data);
 	v->data = NULL;
@@ -51,7 +51,7 @@ struct PsiComplex psi_get_vector(struct PsiVector v, size_t index)
 	return v.data[index];
 }
 
-void psi_set_vector(struct PsiVector *v, size_t index, struct PsiComplex value)
+void psi_set_vector(struct PsiVector* v, size_t index, struct PsiComplex value)
 {
 	assert(index < v->size);
 	v->data[index] = value;
diff --git a/src/psi.c b/src/psi.c
index 6be2831..34940fe 100644
--- a/src/psi.c
+++ b/src/psi.c
@@ -1,6 +1,6 @@
 #include "psi.h"
 
-const char *psi_version(void)
+const char* psi_version(void)
 {
 	return "0.1.0";
 }
diff --git a/src/visualizer/grid.c b/src/visualizer/grid.c
index 5f95efc..32ade1b 100644
--- a/src/visualizer/grid.c
+++ b/src/visualizer/grid.c
@@ -18,7 +18,7 @@ static size_t glyph_bytes(unsigned char lead)
 	return 1;
 }
 
-size_t psi_utf8_count(const char *s)
+size_t psi_utf8_count(const char* s)
 {
 	size_t count = 0;
 	while (*s)
@@ -32,15 +32,14 @@ size_t psi_utf8_count(const char *s)
 
 struct PsiStringBuilder psi_new_string_builder(void)
 {
-	return (struct PsiStringBuilder)
-	{
+	return (struct PsiStringBuilder){
 		NULL,
 		0,
 		0,
 	};
 }
 
-void psi_string_builder_append(struct PsiStringBuilder *sb, const char *s)
+void psi_string_builder_append(struct PsiStringBuilder* sb, const char* s)
 {
 	size_t n = strlen(s);
 	if (sb->len + n + 1 > sb->cap)
@@ -58,7 +57,7 @@ void psi_string_builder_append(struct PsiStringBuilder *sb, const char *s)
 	sb->len += n;
 }
 
-char *psi_string_builder_finish(struct PsiStringBuilder *sb)
+char* psi_string_builder_finish(struct PsiStringBuilder* sb)
 {
 	if (sb->data == NULL)
 		psi_string_builder_append(sb, "");
@@ -74,21 +73,20 @@ struct PsiGlyphRow psi_new_glyph_row(size_t width)
 	for (size_t i = 0; i < width; i++)
 		memcpy(cells[i], " ", 2);
 
-	return (struct PsiGlyphRow)
-	{
+	return (struct PsiGlyphRow){
 		cells,
 		width,
 	};
 }
 
-void psi_free_glyph_row(struct PsiGlyphRow *row)
+void psi_free_glyph_row(struct PsiGlyphRow* row)
 {
 	free(row->cells);
 	row->cells = NULL;
 	row->width = 0;
 }
 
-void psi_glyph_row_set(struct PsiGlyphRow *row, size_t index, const char *glyph)
+void psi_glyph_row_set(struct PsiGlyphRow* row, size_t index, const char* glyph)
 {
 	if (index >= row->width)
 		return;
@@ -98,10 +96,10 @@ void psi_glyph_row_set(struct PsiGlyphRow *row, size_t index, const char *glyph)
 	row->cells[index][n] = '\0';
 }
 
-size_t psi_glyph_row_place(struct PsiGlyphRow *row, size_t start, const char *utf8)
+size_t psi_glyph_row_place(struct PsiGlyphRow* row, size_t start, const char* utf8)
 {
 	size_t count = 0;
-	const char *p = utf8;
+	const char* p = utf8;
 
 	while (*p)
 	{
@@ -120,14 +118,14 @@ size_t psi_glyph_row_place(struct PsiGlyphRow *row, size_t start, const char *ut
 	return count;
 }
 
-void psi_glyph_row_fill_space(struct PsiGlyphRow *row, size_t from, size_t to, const char *glyph)
+void psi_glyph_row_fill_space(struct PsiGlyphRow* row, size_t from, size_t to, const char* glyph)
 {
 	for (size_t i = from; i < to && i < row->width; i++)
 		if (strcmp(row->cells[i], " ") == 0)
 			psi_glyph_row_set(row, i, glyph);
 }
 
-void psi_glyph_row_render(struct PsiGlyphRow row, struct PsiStringBuilder *sb)
+void psi_glyph_row_render(struct PsiGlyphRow row, struct PsiStringBuilder* sb)
 {
 	for (size_t i = 0; i < row.width; i++)
 		psi_string_builder_append(sb, row.cells[i]);
diff --git a/src/visualizer/grid.h b/src/visualizer/grid.h
index f327df8..9353433 100644
--- a/src/visualizer/grid.h
+++ b/src/visualizer/grid.h
@@ -4,14 +4,14 @@
 
 struct PsiStringBuilder
 {
-	char *data;
+	char* data;
 	size_t len;
 	size_t cap;
 };
 
 struct PsiStringBuilder psi_new_string_builder(void);
-void psi_string_builder_append(struct PsiStringBuilder *sb, const char *s);
-char *psi_string_builder_finish(struct PsiStringBuilder *sb);
+void psi_string_builder_append(struct PsiStringBuilder* sb, const char* s);
+char* psi_string_builder_finish(struct PsiStringBuilder* sb);
 
 struct PsiGlyphRow
 {
@@ -20,10 +20,10 @@ struct PsiGlyphRow
 };
 
 struct PsiGlyphRow psi_new_glyph_row(size_t width);
-void psi_free_glyph_row(struct PsiGlyphRow *row);
-void psi_glyph_row_set(struct PsiGlyphRow *row, size_t index, const char *glyph);
-size_t psi_glyph_row_place(struct PsiGlyphRow *row, size_t start, const char *utf8);
-void psi_glyph_row_fill_space(struct PsiGlyphRow *row, size_t from, size_t to, const char *glyph);
-void psi_glyph_row_render(struct PsiGlyphRow row, struct PsiStringBuilder *sb);
+void psi_free_glyph_row(struct PsiGlyphRow* row);
+void psi_glyph_row_set(struct PsiGlyphRow* row, size_t index, const char* glyph);
+size_t psi_glyph_row_place(struct PsiGlyphRow* row, size_t start, const char* utf8);
+void psi_glyph_row_fill_space(struct PsiGlyphRow* row, size_t from, size_t to, const char* glyph);
+void psi_glyph_row_render(struct PsiGlyphRow row, struct PsiStringBuilder* sb);
 
-size_t psi_utf8_count(const char *s);
+size_t psi_utf8_count(const char* s);
diff --git a/src/visualizer/horizontal_cli.c b/src/visualizer/horizontal_cli.c
index e8d0570..56a2e00 100644
--- a/src/visualizer/horizontal_cli.c
+++ b/src/visualizer/horizontal_cli.c
@@ -8,7 +8,7 @@
 
 #include "visualizer/grid.h"
 
-static void append_repeat(struct PsiStringBuilder *sb, const char *glyph, size_t n)
+static void append_repeat(struct PsiStringBuilder* sb, const char* glyph, size_t n)
 {
 	for (size_t i = 0; i < n; i++)
 		psi_string_builder_append(sb, glyph);
@@ -34,14 +34,12 @@ static bool is_single_target(enum PsiGateKind kind)
 		case PSI_GATE_P:
 		case PSI_GATE_U1:
 		case PSI_GATE_U2:
-		case PSI_GATE_U3:
-			return true;
-		default:
-			return false;
+		case PSI_GATE_U3: return true;
+		default: return false;
 	}
 }
 
-static void single_label(struct PsiGateOp op, char *out, size_t cap)
+static void single_label(struct PsiGateOp op, char* out, size_t cap)
 {
 	switch (op.kind)
 	{
@@ -66,7 +64,7 @@ static void single_label(struct PsiGateOp op, char *out, size_t cap)
 	}
 }
 
-static void controlled_label(struct PsiGateOp op, char *out, size_t cap)
+static void controlled_label(struct PsiGateOp op, char* out, size_t cap)
 {
 	switch (op.kind)
 	{
@@ -85,20 +83,18 @@ static bool is_param_controlled(enum PsiGateKind kind)
 		case PSI_GATE_CRX:
 		case PSI_GATE_CRY:
 		case PSI_GATE_CRZ:
-		case PSI_GATE_CP:
-			return true;
-		default:
-			return false;
+		case PSI_GATE_CP: return true;
+		default: return false;
 	}
 }
 
-char *psi_render_circuit_horizontal(const struct PsiQuantumCircuit *circuit)
+char* psi_render_circuit_horizontal(const struct PsiQuantumCircuit* circuit)
 {
 	size_t nq = circuit->num_qubits;
 	size_t nc = circuit->num_classical;
 
-	struct PsiStringBuilder *q = nq > 0 ? malloc(nq * sizeof(struct PsiStringBuilder)) : NULL;
-	struct PsiStringBuilder *c = nc > 0 ? malloc(nc * sizeof(struct PsiStringBuilder)) : NULL;
+	struct PsiStringBuilder* q = nq > 0 ? malloc(nq * sizeof(struct PsiStringBuilder)) : NULL;
+	struct PsiStringBuilder* c = nc > 0 ? malloc(nc * sizeof(struct PsiStringBuilder)) : NULL;
 	assert(q != NULL || nq == 0);
 	assert(c != NULL || nc == 0);
 
@@ -152,7 +148,7 @@ char *psi_render_circuit_horizontal(const struct PsiQuantumCircuit *circuit)
 		struct PsiGateOp op = circuit->operations[oi];
 
 		size_t tc;
-		const size_t *targets = psi_gate_op_quantum_targets(&op, &tc);
+		const size_t* targets = psi_gate_op_quantum_targets(&op, &tc);
 		size_t min_q = targets[0];
 		size_t max_q = targets[0];
 		for (size_t k = 1; k < tc; k++)
@@ -216,40 +212,52 @@ char *psi_render_circuit_horizontal(const struct PsiQuantumCircuit *circuit)
 				append_repeat(&c[i], "═", w + 2);
 			append_repeat(&gap, " ", w + 2);
 		}
-		else if (op.kind == PSI_GATE_CNOT || op.kind == PSI_GATE_CZ || op.kind == PSI_GATE_SWAP
-			|| op.kind == PSI_GATE_CCNOT || op.kind == PSI_GATE_CSWAP)
+		else if (op.kind == PSI_GATE_CNOT || op.kind == PSI_GATE_CZ || op.kind == PSI_GATE_SWAP ||
+		         op.kind == PSI_GATE_CCNOT || op.kind == PSI_GATE_CSWAP)
 		{
 			for (size_t i = 0; i < nq; i++)
 			{
-				const char *seg = "─────";
+				const char* seg = "─────";
 
 				switch (op.kind)
 				{
 					case PSI_GATE_CNOT:
-						if (i == op.qubits[0]) seg = "──●──";
-						else if (i == op.qubits[1]) seg = "──⊕──";
-						else if (i > min_q && i < max_q) seg = "──│──";
+						if (i == op.qubits[0])
+							seg = "──●──";
+						else if (i == op.qubits[1])
+							seg = "──⊕──";
+						else if (i > min_q && i < max_q)
+							seg = "──│──";
 						break;
 					case PSI_GATE_CZ:
-						if (i == op.qubits[0] || i == op.qubits[1]) seg = "──●──";
-						else if (i > min_q && i < max_q) seg = "──│──";
+						if (i == op.qubits[0] || i == op.qubits[1])
+							seg = "──●──";
+						else if (i > min_q && i < max_q)
+							seg = "──│──";
 						break;
 					case PSI_GATE_SWAP:
-						if (i == op.qubits[0] || i == op.qubits[1]) seg = "──╳──";
-						else if (i > min_q && i < max_q) seg = "──│──";
+						if (i == op.qubits[0] || i == op.qubits[1])
+							seg = "──╳──";
+						else if (i > min_q && i < max_q)
+							seg = "──│──";
 						break;
 					case PSI_GATE_CCNOT:
-						if (i == op.qubits[0] || i == op.qubits[1]) seg = "──●──";
-						else if (i == op.qubits[2]) seg = "──⊕──";
-						else if (i > min_q && i < max_q) seg = "──│──";
+						if (i == op.qubits[0] || i == op.qubits[1])
+							seg = "──●──";
+						else if (i == op.qubits[2])
+							seg = "──⊕──";
+						else if (i > min_q && i < max_q)
+							seg = "──│──";
 						break;
 					case PSI_GATE_CSWAP:
-						if (i == op.qubits[0]) seg = "──●──";
-						else if (i == op.qubits[1] || i == op.qubits[2]) seg = "──╳──";
-						else if (i > min_q && i < max_q) seg = "──│──";
-						break;
-					default:
+						if (i == op.qubits[0])
+							seg = "──●──";
+						else if (i == op.qubits[1] || i == op.qubits[2])
+							seg = "──╳──";
+						else if (i > min_q && i < max_q)
+							seg = "──│──";
 						break;
+					default: break;
 				}
 
 				psi_string_builder_append(&q[i], seg);
diff --git a/src/visualizer/vertical_cli.c b/src/visualizer/vertical_cli.c
index 09d4456..4f21e49 100644
--- a/src/visualizer/vertical_cli.c
+++ b/src/visualizer/vertical_cli.c
@@ -5,7 +5,7 @@
 
 #include "visualizer/grid.h"
 
-static void gate_label(struct PsiGateOp op, char *out, size_t cap)
+static void gate_label(struct PsiGateOp op, char* out, size_t cap)
 {
 	switch (op.kind)
 	{
@@ -40,7 +40,7 @@ static void gate_label(struct PsiGateOp op, char *out, size_t cap)
 	}
 }
 
-static size_t calculate_col_width(const struct PsiQuantumCircuit *circuit)
+static size_t calculate_col_width(const struct PsiQuantumCircuit* circuit)
 {
 	size_t max_label = 3;
 
@@ -80,10 +80,8 @@ static bool is_single_target(enum PsiGateKind kind)
 		case PSI_GATE_P:
 		case PSI_GATE_U1:
 		case PSI_GATE_U2:
-		case PSI_GATE_U3:
-			return true;
-		default:
-			return false;
+		case PSI_GATE_U3: return true;
+		default: return false;
 	}
 }
 
@@ -94,14 +92,12 @@ static bool is_param_controlled(enum PsiGateKind kind)
 		case PSI_GATE_CRX:
 		case PSI_GATE_CRY:
 		case PSI_GATE_CRZ:
-		case PSI_GATE_CP:
-			return true;
-		default:
-			return false;
+		case PSI_GATE_CP: return true;
+		default: return false;
 	}
 }
 
-char *psi_render_circuit_vertical(const struct PsiQuantumCircuit *circuit)
+char* psi_render_circuit_vertical(const struct PsiQuantumCircuit* circuit)
 {
 	struct PsiStringBuilder sb = psi_new_string_builder();
 
@@ -154,7 +150,7 @@ char *psi_render_circuit_vertical(const struct PsiQuantumCircuit *circuit)
 		struct PsiGateOp op = circuit->operations[op_index];
 
 		size_t target_count;
-		const size_t *targets = psi_gate_op_quantum_targets(&op, &target_count);
+		const size_t* targets = psi_gate_op_quantum_targets(&op, &target_count);
 		size_t min_q = targets[0];
 		size_t max_q = targets[0];
 		for (size_t i = 1; i < target_count; i++)
@@ -212,8 +208,8 @@ char *psi_render_circuit_vertical(const struct PsiQuantumCircuit *circuit)
 		else if (op.kind == PSI_GATE_CCNOT || op.kind == PSI_GATE_CSWAP)
 		{
 			bool is_cswap = op.kind == PSI_GATE_CSWAP;
-			const char *sym_c = "●";
-			const char *sym_t = is_cswap ? "╳" : "⊕";
+			const char* sym_c = "●";
+			const char* sym_t = is_cswap ? "╳" : "⊕";
 			size_t c1 = op.qubits[0];
 			size_t c2 = op.qubits[1];
 			size_t t = op.qubits[2];
@@ -267,7 +263,8 @@ char *psi_render_circuit_vertical(const struct PsiQuantumCircuit *circuit)
 						psi_glyph_row_set(&row, center, "□");
 				}
 
-				psi_glyph_row_fill_space(&row, min_q * stride + half + 1, max_q * stride + half, "─");
+				psi_glyph_row_fill_space(&row, min_q * stride + half + 1, max_q * stride + half,
+				                         "─");
 			}
 
 			for (size_t i = 0; i < nc; i++)
@@ -275,8 +272,8 @@ char *psi_render_circuit_vertical(const struct PsiQuantumCircuit *circuit)
 		}
 		else
 		{
-			const char *sym1 = "●";
-			const char *sym2 = "⊕";
+			const char* sym1 = "●";
+			const char* sym2 = "⊕";
 			if (op.kind == PSI_GATE_CZ)
 				sym2 = "●";
 			else if (op.kind == PSI_GATE_SWAP)
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);