diff --git a/circuit.lisp b/circuit.lisp
index 5ef2d90..348c13c 100644
--- a/circuit.lisp
+++ b/circuit.lisp
@@ -4,13 +4,13 @@
 (defparameter *circuit-measure-places* 5
   "The number of places to use with the :measure instruction in circuits.")
 
-(defun circuit-arbitrary-gate (state operator gate target &optional control)
+(defun circuit-arbitrary-gate (size operator gate target &optional control)
   (declare (ignorable operator))
   (if control
-      (napply-controlled-gate state gate target control)
-      (napply-gate state gate target)))
+      (make-controlled-operator size gate target control)
+      (make-operator size gate target)))
 
-(defun circuit-matrix-operation (state operator target &optional control)
+(defun circuit-matrix-operation (size operator target &optional control)
   (let ((matrix (case operator
                   (:* +identity-2x2+)
                   (:x +pauli-x-gate+)
@@ -22,16 +22,13 @@
                   (:cnot +pauli-x-gate+)
                   (:cz +pauli-z-gate+))))
     (if control
-        (napply-controlled-gate state matrix target control)
-        (napply-gate state matrix target))))
+        (make-controlled-operator size matrix target control)
+        (make-operator size matrix target))))
 
-(defun circuit-n-controlled-gate (state operator target controls)
-  (let ((func (case operator
-                (:ntoff 'make-n-toffoli-operator)
-                (:ncz 'make-n-controlled-z-operator))))
-    (replace-state state (*mv (funcall func (state-bits state)
-                                       target controls)
-                              state))))
+(defun circuit-n-controlled-gate (size operator target controls)
+  (case operator
+    (:ntoff (make-n-toffoli-operator size target controls))
+    (:ncz (make-n-controlled-z-operator size target controls))))
 
 (defparameter *circuit-operators*
   ;; Operator, # args, more args?, has output?, function
@@ -47,36 +44,72 @@
     (:t 1 nil nil circuit-matrix-operation)
     (:cnot 2 nil nil circuit-matrix-operation)
     (:cz 2 nil nil circuit-matrix-operation)
-    (:measure 2 nil t ,(lambda (state operator &rest args)
-                         (declare (ignorable operator))
-                         (nmeasure state (car args)
-                                   :places *circuit-measure-places*)))
     (:ntoff 2 nil nil circuit-n-controlled-gate)
-    (:ncz 2 nil nil circuit-n-controlled-gate)))
+    (:ncz 2 nil nil circuit-n-controlled-gate)
+    (:measure 2 nil t ,(lambda (size operator bit target)
+                         (declare (ignorable size operator target))
+                         (lambda (state)
+                           (nmeasure state bit
+                                     :places *circuit-measure-places*))))))
 
-(defun make-circuit ()
-  "Create a new blank circuit."
-  '(:circuit))
+(defun make-circuit (bits)
+  "Create a new blank circuit of BITS bits."
+  `(:circuit ,bits))
+
+(defun circuitp (obj)
+  "Return non-nil if OBJ is a quantum circuit."
+  (and (eq (first obj) :circuit)
+       (integerp (second obj))))
 
 (defun add-to-circuit (circuit operator &rest args)
   "Add OPERATOR to CIRCUIT."
-  (let ((entry (assoc operator *circuit-operators*)))
+  (let ((bits (second circuit))
+        (entry (assoc operator *circuit-operators*)))
     (unless entry
       (error "Unknown circuit operator: ~s" operator))
-    (destructuring-bind (name arg-count has-rest &rest r) entry
-      (declare (ignorable name r))
+    (destructuring-bind (name arg-count has-rest has-output func) entry
+      (declare (ignorable name))
       (unless (or (and has-rest (>= (length args) arg-count))
                   (= (length args) arg-count))
         (error "Operator ~s expects ~@[~*exactly ~]~s arg~:p, got ~s" operator
                (not has-rest) arg-count (length args)))
-      (nconc circuit (list (cons operator args))))
+      (let ((action (apply func bits operator args))
+            (prev-action (car (last circuit))))
+        (cond
+          ((and (matrixp action) (matrixp prev-action))
+           (setf (car (last circuit)) (*mm action prev-action)))
+          ((functionp action)
+           (setf (cdr (last circuit)) (list (cons action
+                                                  (and has-output
+                                                       (car (last args)))))))
+          (t (setf (cdr (last circuit)) (list action))))))
     circuit))
 
-(defmacro with-build-circuit (&body body)
+(defun build-circuit-add-action (bits actions args)
+  "Add the action denoted by ARGS to the front of ACTIONS. The action is assumed
+to operate on a state of BITS bits."
+  (let ((entry (assoc (car args) *circuit-operators*)))
+    (destructuring-bind (name arg-count has-rest has-output func) entry
+      (declare (ignorable arg-count has-rest))
+      (declare (ignorable name))
+      (let ((action (apply func bits (car args) (cdr args)))
+            (prev-action (car actions)))
+        (cond
+          ((and (matrixp action) (matrixp prev-action))
+           (setf (car actions) (*mm action prev-action)))
+          ((functionp action)
+           (push (cons action (and has-output (car (last args))))
+                 actions))
+          (t (push action actions)))
+        actions))))
+
+(defmacro with-build-circuit (bits &body body)
   "Create a circuit using a simple DSL. BODY can be any valid Lisp forms, in
 addition to function calls to functions named in `*circuit-operators*'."
-  (let ((circuit-var (gensym)))
-    `(let ((,circuit-var))
+  (let ((circuit-var (gensym))
+        (bits-var (gensym)))
+    `(let ((,circuit-var)
+           (,bits-var ,bits))
        (macrolet
            (,@(mapcar (lambda (oper)
                         (let ((arg (gensym))
@@ -91,44 +124,37 @@ addition to function calls to functions named in `*circuit-operators*'."
                             (declare (ignorable ,@arg-list
                                                 ,@(when (third oper)
                                                     (list whole-arg))))
-                            `(push (list ,@,arg) ,',circuit-var))))
+                            `(setq ,',circuit-var
+                                   (build-circuit-add-action
+                                    ,',bits-var ,',circuit-var (list ,@,arg))))))
                 *circuit-operators*))
          ,@body)
-       (cons :circuit (nreverse ,circuit-var)))))
+       (apply 'list :circuit ,bits-var (nreverse ,circuit-var)))))
 
-(defun apply-circuit-operator-to-state (state operator args)
-  "Apply the circuit operator OPERATOR to STATE by calling its function with
-ARGS."
-  (destructuring-bind (&optional name arg-count has-rest has-output function)
-      (assoc operator *circuit-operators*)
-    (declare (ignorable name arg-count has-rest))
-    (assert function ()
-            "Invalid circuit operator: ~s" operator)
-    (let ((output (apply function state operator args)))
-      (when has-output
-        (cons (car (last args)) output)))))
-
-(defun run-circuit (circuit &key bits uniform coefficients probabilities (places 5))
-  "Run the circuit CIRCUIT and return the final state. The initial state can be
-specified in one of three ways:
- - BITS: the number of qbits
- - COEFFICIENTS: the initial coefficients
- - PROBABILITES: the initial probabilities"
-  (assert (= 1 (count-if 'identity (list bits coefficients probabilities uniform)))
-          ()
-          "Exactly one of BITS, UNIFORM, COEFFICIENTS, and PROBABILITIES can ~
-be present")
-  (let ((state (cond
-                 (bits (make-zero-state bits))
-                 (uniform (make-uniform-normal-state uniform))
-                 (coefficients (coerce coefficients 'vector))
-                 (probabilities (make-normal-state probabilities)))))
+(defun run-circuit (circuit &key (state :zero) (places 5))
+  "Run the circuit CIRCUIT and return the final state. The initial STATE can be
+specified in one of four ways:
+ - :ZERO: a zero state is used -> #(1 0 0 0), etc.
+ - :UNIFORM: a uniform state -> #(1/sqrt(2) 1/sqrt(2)), etc
+ - a vector: the initial coefficients"
+  (let ((real-state (cond
+                 ((eq state :zero)
+                  (make-zero-state (second circuit)))
+                 ((eq state :uniform)
+                  (make-uniform-normal-state (second circuit)))
+                 (t (coerce state 'vector)))))
     (values
-     state
+     real-state
      (loop with *circuit-measure-places* = places
-           for element in (cdr circuit)
-           for name = (car element)
-           for args = (cdr element)
-           for result = (apply-circuit-operator-to-state state name args)
-           when result
-             collect result))))
+           with outputs = ()
+           for action in (cddr circuit) do
+             (cond
+               ((matrixp action)
+                (replace-state real-state (*mv action real-state)))
+               ((consp action)
+                (destructuring-bind (func output) action
+                  (let ((res (funcall func real-state)))
+                    (when output
+                      (push res outputs)
+                      (push output outputs)))))
+               (t (error "Unknown circuit element type: ~s" action)))))))
diff --git a/examples/bell.lisp b/examples/bell.lisp
index 34efbc9..5112363 100644
--- a/examples/bell.lisp
+++ b/examples/bell.lisp
@@ -1,29 +1,29 @@
 ;;;; Example Bell state creation circuits
 (in-package :cl-quantum/examples)
 
-(defun make-bell-phi+-circuit ()
+(defun make-bell-phi+-circuit (bits)
   "Return a circuit that makes a Bell phi+ state. That is |00>+|11>."
-  (with-build-circuit
+  (with-build-circuit bits
     (:h 0)
     (:cnot 1 0)))
 
-(defun make-bell-phi--circuit ()
+(defun make-bell-phi--circuit (bits)
   "Return a circuit that makes a Bell phi- state. That is |00>-|11>."
-  (with-build-circuit
+  (with-build-circuit bits
     (:h 0)
     (:cnot 1 0)
     (:z 1)))
 
-(defun make-bell-psi+-circuit ()
+(defun make-bell-psi+-circuit (bits)
   "Return a circuit that makes a Bell psi+ state. That is |01>+|10>."
-  (with-build-circuit
+  (with-build-circuit bits
     (:x 1)
     (:h 0)
     (:cnot 1 0)))
 
-(defun make-bell-psi--circuit ()
+(defun make-bell-psi--circuit (bits)
   "Return a circuit that makes a Bell psi- state. That is |01>-|10>."
-  (with-build-circuit
+  (with-build-circuit bits
     (:x 1)
     (:h 0)
     (:cnot 1 0)
diff --git a/examples/grover.lisp b/examples/grover.lisp
index e185662..1213064 100644
--- a/examples/grover.lisp
+++ b/examples/grover.lisp
@@ -12,7 +12,7 @@ bits and finds when the state is equal to TARGET."
   (assert (> (ash 1 bits) target)
           (bits target)
           "Target bit of ~s out of range for state with ~s bits." target bits)
-  (with-build-circuit
+  (with-build-circuit bits
     ;; Setup
     (loop for i below bits do (:h i))
 
diff --git a/math.lisp b/math.lisp
index 945423a..cd2268a 100644
--- a/math.lisp
+++ b/math.lisp
@@ -42,6 +42,11 @@ VALUE, the ROW, and the COLUMN."
    (domatrix ((elem row col) matrix new-mat)
      (setf (aref new-mat row col) (funcall function elem row col)))))
 
+(defun matrixp (obj)
+  "Return non-nil if OBJ is a matrix."
+  (and (arrayp obj)
+       (= (array-rank obj) 2)))
+
 ;; Matrix subroutines
 (defun minor (mat i j)
   "Find the minor of MAT for I and J."
diff --git a/package.lisp b/package.lisp
index 714220e..604a768 100644
--- a/package.lisp
+++ b/package.lisp
@@ -4,6 +4,7 @@ quantum operators.")
   (:use :cl)
   (:export #:domatrix
            #:mapmatrix
+           #:matrixp
            #:minor
            #:cofactor
            #:det