diff --git a/pyproject.toml b/pyproject.toml
index e731b2a..f31b966 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -1,6 +1,6 @@
 [project]
 name = "grovers-visualizer"
-version = "0.4.3"
+version = "0.5.0"
 description = "A tiny Python package that steps through Grover’s Search algorithm."
 readme = "README.md"
 requires-python = ">=3.10"
diff --git a/src/grovers_visualizer/circuit.py b/src/grovers_visualizer/circuit.py
index 352c945..e93d831 100644
--- a/src/grovers_visualizer/circuit.py
+++ b/src/grovers_visualizer/circuit.py
@@ -11,6 +11,13 @@ def oracle(qc: QuantumCircuit, target_state: QubitState) -> None:
     encode_target_state(qc, target_state)  # Undo
 
 
+def oracle_circuit(target: QubitState) -> QuantumCircuit:
+    n = len(target)
+    qc = QuantumCircuit(n)
+    oracle(qc, target)
+    return qc
+
+
 def diffusion(qc: QuantumCircuit, n: int) -> None:
     """Apply the Grovers diffusion operator."""
     qc.h(range(n))
@@ -20,6 +27,12 @@ def diffusion(qc: QuantumCircuit, n: int) -> None:
     qc.h(range(n))
 
 
+def diffusion_circuit(n: int) -> QuantumCircuit:
+    qc = QuantumCircuit(n)
+    diffusion(qc, n)
+    return qc
+
+
 def encode_target_state(qc: QuantumCircuit, target_state: QubitState) -> None:
     """Apply X gates to qubits where the target state bit is '0'."""
     for i, bit in enumerate(reversed(target_state)):
diff --git a/src/grovers_visualizer/simulation.py b/src/grovers_visualizer/simulation.py
index 0bdfc15..f3ca0bc 100644
--- a/src/grovers_visualizer/simulation.py
+++ b/src/grovers_visualizer/simulation.py
@@ -2,30 +2,32 @@ from collections.abc import Iterator
 from itertools import count
 
 from qiskit import QuantumCircuit
+from qiskit.qasm2.parse import Operator
 from qiskit.quantum_info import Statevector
 
-from grovers_visualizer.circuit import diffusion, oracle
+from grovers_visualizer.circuit import diffusion_circuit, oracle_circuit
 from grovers_visualizer.state import QubitState
 
 
 def grover_evolver(target: QubitState, max_iterations: int = 0) -> Iterator[tuple[int, Statevector]]:
     """Yields (iteration, statevector) pairs.
 
-    iteration=0 is the uniform-Hadamard initialization. If
-    max_iterations > 0, stop after that many iterations. If
-    max_iterations == 0, run indefinitely (until the consumer breaks).
+    - iteration=0 is the uniform-Hadamard initialization
+    - max_iterations > 0, stop after that many iterations
+    - max_iterations == 0, run indefinitely (until the consumer breaks)
     """
     n_qubits = len(target)
     qc = QuantumCircuit(n_qubits)
     qc.h(range(n_qubits))
 
     # initial statevector
-    yield 0, Statevector.from_instruction(qc)
+    sv = Statevector.from_instruction(qc)
+    yield 0, sv
 
-    # pick an iterator for subsequent steps
-    iter = range(1, max_iterations + 1) if max_iterations > 0 else count(1)
+    oracle_op = Operator(oracle_circuit(target))
+    diffusion_op = Operator(diffusion_circuit(n_qubits))
 
-    for i in iter:
-        oracle(qc, target)
-        diffusion(qc, n_qubits)
-        yield i, Statevector.from_instruction(qc)
+    iters = range(1, max_iterations + 1) if max_iterations > 0 else count(1)
+    for i in iters:
+        sv = sv.evolve(oracle_op).evolve(diffusion_op)
+        yield i, sv