diff --git a/CHANGELOG b/CHANGELOG
index af238b94b..e482e892f 100644
--- a/CHANGELOG
+++ b/CHANGELOG
@@ -1,5 +1,17 @@
 # CHANGELOG
 
+## [0.9.12.3] - 2024-06-03
+
+### Added
+* Deterministic Clifford compilation and native gate count statistics for `CliffordRBDesign` (#314, #315, #443)
+
+
+### Fixed
+* Truncation bugfix in `BenchmarkingDesign` objects with "paired" lists to `circuit_list` attribute (#408, #443)
+* Fixes and efficiency improvements for various linear algebra calls (#432)
+* `densitymx_slow` evotype hotfix (#438, #439)
+
+
 ## [0.9.12.2] - 2024-04-16
 
 ### Added
diff --git a/pygsti/algorithms/randomcircuit.py b/pygsti/algorithms/randomcircuit.py
index c45f1d76c..ee7140e01 100644
--- a/pygsti/algorithms/randomcircuit.py
+++ b/pygsti/algorithms/randomcircuit.py
@@ -2137,9 +2137,129 @@ def create_direct_rb_circuit(pspec, clifford_compilations, length, qubit_labels=
 
 #     return experiment_dict
 
+def _sample_clifford_circuit(pspec, clifford_compilations, qubit_labels, citerations,
+        compilerargs, exact_compilation_key, srep_cache, rand_state):
+    """Helper function to compile a random Clifford circuit.
+
+    Parameters
+    ----------
+    pspec : QubitProcessorSpec
+        The QubitProcessorSpec for the device that the circuit is being sampled for, which defines the
+        "native" gate-set and the connectivity of the device. The returned CRB circuit will be over
+        the gates in `pspec`, and will respect the connectivity encoded by `pspec`.
+
+    clifford_compilations : dict
+        A dictionary with at least the potential keys `'absolute'` and `'paulieq'` and corresponding
+        :class:`CompilationRules` values.  These compilation rules specify how to compile the
+        "native" gates of `pspec` into Clifford gates. Additional :class:`CompilationRules` can be
+        provided, particularly for use with `exact_compilation_key`.
+
+    qubit_labels : list
+        A list of the qubits that the RB circuit is to be sampled for.
+
+    citerations : int
+        Some of the Clifford compilation algorithms in pyGSTi (including the default algorithm) are
+        randomized, and the lowest-cost circuit is chosen from all the circuit generated in the
+        iterations of the algorithm. This is the number of iterations used. The time required to
+        generate a CRB circuit is linear in `citerations` * (`length` + 2). Lower-depth / lower 2-qubit
+        gate count compilations of the Cliffords are important in order to successfully implement
+        CRB on more qubits.
+
+    compilerargs : list
+        A list of arguments that are handed to compile_clifford() function, which includes all the
+        optional arguments of compile_clifford() *after* the `iterations` option (set by `citerations`).
+        In order, this list should be values for:
+        
+        algorithm : str. A string that specifies the compilation algorithm. The default in
+        compile_clifford() will always be whatever we consider to be the 'best' all-round
+        algorithm
+        
+        aargs : list. A list of optional arguments for the particular compilation algorithm.
+        
+        costfunction : 'str' or function. The cost-function from which the "best" compilation
+        for a Clifford is chosen from all `citerations` compilations. The default costs a
+        circuit as 10x the num. of 2-qubit gates in the circuit + 1x the depth of the circuit.
+        
+        prefixpaulis : bool. Whether to prefix or append the Paulis on each Clifford.
+        
+        paulirandomize : bool. Whether to follow each layer in the Clifford circuit with a
+        random Pauli on each qubit (compiled into native gates). I.e., if this is True the
+        native gates are Pauli-randomized. When True, this prevents any coherent errors adding
+        (on average) inside the layers of each compiled Clifford, at the cost of increased
+        circuit depth. Defaults to False.
+        
+        For more information on these options, see the `:func:compile_clifford()` docstring.
+    
+    exact_compilation_key: str, optional
+        The key into `clifford_compilations` to use for exact deterministic complation of Cliffords.
+        The underlying :class:`CompilationRules` object must provide compilations for all possible
+        n-qubit Cliffords that will be generated. This also requires the pspec is able to generate the
+        symplectic representations for all n-qubit Cliffords in :meth:`compute_clifford_symplectic_reps`.
+        This is currently generally intended for use out-of-the-box with 1-qubit Clifford RB;
+        however, larger number of qubits can be used so long as the user specifies the processor spec and
+        compilation rules properly.
+
+    srep_cache: dict
+        Keys are gate labels and values are precomputed symplectic representations.
+    
+    rand_state: np.random.RandomState
+        A RandomState to use for RNG
+
+    Returns
+    -------
+    clifford_circuit : Circuit
+        The compiled Clifford circuit
+    
+    s:
+        The symplectic matrix of the Clifford
+    
+    p:
+        The symplectic phase vector of the Clifford
+    """
+    # Find the labels of the qubits to create the circuit for.
+    if qubit_labels is not None: qubits = qubit_labels[:]  # copy this list
+    else: qubits = pspec.qubit_labels[:]  # copy this list
+    # The number of qubits the circuit is over.
+    n = len(qubits)
+
+    if exact_compilation_key is not None:
+        # Deterministic compilation based on a provided clifford compilation
+        assert exact_compilation_key in clifford_compilations, \
+                f"{exact_compilation_key} not provided in `clifford_compilations`"
+        
+        # Pick clifford
+        cidx = rand_state.randint(_symp.compute_num_cliffords(n))
+        lbl = _lbl.Label(f'C{cidx}', qubits)
+        
+        # Try to do deterministic compilation
+        try:
+            circuit = clifford_compilations[exact_compilation_key].retrieve_compilation_of(lbl)
+        except AssertionError:
+            raise ValueError(
+                f"Failed to compile n-qubit Clifford 'C{cidx}'. Ensure this is provided in the " + \
+                "compilation rules, or use a compilation algorithm to synthesize it by not " + \
+                "specifying `exact_compilation_key`."
+            )
+    
+        # compute the symplectic rep of the chosen clifford
+        # TODO: Note that this is inefficient. For speed, we could implement the pair to
+        # _symp.compute_symplectic_matrix and just calculate s and p directly
+        s, p = _symp.symplectic_rep_of_clifford_circuit(circuit, srep_cache)
+    else:
+        # Random compilation
+        s, p = _symp.random_clifford(n, rand_state=rand_state)
+        circuit = _cmpl.compile_clifford(s, p, pspec,
+                                        clifford_compilations.get('absolute', None),
+                                        clifford_compilations.get('paulieq', None),
+                                        qubit_labels=qubit_labels, iterations=citerations, *compilerargs,
+                                            rand_state=rand_state)
+    
+    return circuit, s, p
+
 
 def create_clifford_rb_circuit(pspec, clifford_compilations, length, qubit_labels=None, randomizeout=False,
-                               citerations=20, compilerargs=None, interleaved_circuit=None, seed=None):
+                               citerations=20, compilerargs=None, interleaved_circuit=None, seed=None,
+                               return_native_gate_counts=False, exact_compilation_key=None):
     """
     Generates a "Clifford randomized benchmarking" (CRB) circuit.
 
@@ -2165,9 +2285,10 @@ def create_clifford_rb_circuit(pspec, clifford_compilations, length, qubit_label
         the gates in `pspec`, and will respect the connectivity encoded by `pspec`.
 
     clifford_compilations : dict
-        A dictionary with the potential keys `'absolute'` and `'paulieq'` and corresponding
+        A dictionary with at least the potential keys `'absolute'` and `'paulieq'` and corresponding
         :class:`CompilationRules` values.  These compilation rules specify how to compile the
-        "native" gates of `pspec` into Clifford gates.
+        "native" gates of `pspec` into Clifford gates. Additional :class:`CompilationRules` can be
+        provided, particularly for use with `exact_compilation_key`.
 
     length : int
         The "CRB length" of the circuit -- an integer >= 0 --  which is the number of Cliffords in the
@@ -2223,6 +2344,18 @@ def create_clifford_rb_circuit(pspec, clifford_compilations, length, qubit_label
     seed : int, optional
         A seed to initialize the random number generator used for creating random clifford
         circuits.
+    
+    return_native_gate_counts: bool, optional
+        Whether to return the number of native gates in the first `length`+1 compiled Cliffords
+    
+    exact_compilation_key: str, optional
+        The key into `clifford_compilations` to use for exact deterministic complation of Cliffords.
+        The underlying :class:`CompilationRules` object must provide compilations for all possible
+        n-qubit Cliffords that will be generated. This also requires the pspec is able to generate the
+        symplectic representations for all n-qubit Cliffords in :meth:`compute_clifford_symplectic_reps`.
+        This is currently generally intended for use out-of-the-box with 1-qubit Clifford RB;
+        however, larger number of qubits can be used so long as the user specifies the processor spec and
+        compilation rules properly.
 
     Returns
     -------
@@ -2236,6 +2369,10 @@ def create_clifford_rb_circuit(pspec, clifford_compilations, length, qubit_label
         `qubit_labels`, if `qubit_labels` is not None; the ith element of `pspec.qubit_labels`, otherwise.
         In both cases, the ith element of the tuple corresponds to the error-free outcome for the
         qubit on the ith wire of the output circuit.
+    
+    native_gate_counts: dict
+        Total number of native gates, native 2q gates, and native circuit size in the
+        first `length`+1 compiled Cliffords. Only returned when `return_num_native_gates` is True
     """
     if compilerargs is None:
         compilerargs = []
@@ -2245,6 +2382,12 @@ def create_clifford_rb_circuit(pspec, clifford_compilations, length, qubit_label
     # The number of qubits the circuit is over.
     n = len(qubits)
 
+    srep_cache = {}
+    if exact_compilation_key is not None:
+        # Precompute some of the symplectic reps if we are doing exact compilation
+        srep_cache = _symp.compute_internal_gate_symplectic_representations()
+        srep_cache.update(pspec.compute_clifford_symplectic_reps())
+
     rand_state = _np.random.RandomState(seed)  # OK if seed is None
 
     # Initialize the identity circuit rep.
@@ -2255,14 +2398,17 @@ def create_clifford_rb_circuit(pspec, clifford_compilations, length, qubit_label
 
     # Sample length+1 uniformly random Cliffords (we want a circuit of length+2 Cliffords, in total), compile
     # them, and append them to the current circuit.
-    for i in range(0, length + 1):
+    num_native_gates = 0
+    num_native_2q_gates = 0
+    native_size = 0
+    for _ in range(0, length + 1):
+        # Perform sampling
+        circuit, s, p = _sample_clifford_circuit(pspec, clifford_compilations, qubit_labels, citerations,
+                                 compilerargs, exact_compilation_key, srep_cache, rand_state)
+        num_native_gates += circuit.num_gates
+        num_native_2q_gates += circuit.num_nq_gates(2)
+        native_size += circuit.size
 
-        s, p = _symp.random_clifford(n, rand_state=rand_state)
-        circuit = _cmpl.compile_clifford(s, p, pspec,
-                                         clifford_compilations.get('absolute', None),
-                                         clifford_compilations.get('paulieq', None),
-                                         qubit_labels=qubit_labels, iterations=citerations, *compilerargs,
-                                         rand_state=rand_state)
         # Keeps track of the current composite Clifford
         s_composite, p_composite = _symp.compose_cliffords(s_composite, p_composite, s, p)
         full_circuit.append_circuit_inplace(circuit)
@@ -2306,6 +2452,15 @@ def create_clifford_rb_circuit(pspec, clifford_compilations, length, qubit_label
     idealout = tuple(idealout)
 
     full_circuit.done_editing()
+
+    native_gate_counts = {
+        "native_gate_count": num_native_gates,
+        "native_2q_gate_count": num_native_2q_gates,
+        "native_size": native_size
+    }
+
+    if return_native_gate_counts:
+        return full_circuit, idealout, native_gate_counts
     return full_circuit, idealout
 
 
diff --git a/pygsti/circuits/circuit.py b/pygsti/circuits/circuit.py
index 8e629031a..25cbdd923 100644
--- a/pygsti/circuits/circuit.py
+++ b/pygsti/circuits/circuit.py
@@ -3405,6 +3405,30 @@ def two_q_gate_count(self):
         """
         return self.num_nq_gates(2)
 
+    @property
+    def num_gates(self):
+        """
+        The number of gates in the circuit.
+
+        Returns
+        -------
+        int
+        """
+        if self._static:
+            def cnt(lbl):  # obj a Label, perhaps compound
+                if lbl.is_simple():  # a simple label
+                    return 1 if (lbl.sslbls is not None) else 0
+                else:
+                    return sum([cnt(sublbl) for sublbl in lbl.components])
+        else:
+            def cnt(obj):  # obj is either a simple label or a list
+                if isinstance(obj, _Label):  # all Labels are simple labels
+                    return 1 if (obj.sslbls is not None) else 0
+                else:
+                    return sum([cnt(sub) for sub in obj])
+
+        return sum([cnt(layer_lbl) for layer_lbl in self._labels])
+
     def num_nq_gates(self, nq):
         """
         The number of `nq`-qubit gates in the circuit.
diff --git a/pygsti/circuits/circuitlist.py b/pygsti/circuits/circuitlist.py
index 3c9345269..666802382 100644
--- a/pygsti/circuits/circuitlist.py
+++ b/pygsti/circuits/circuitlist.py
@@ -158,11 +158,7 @@ def truncate(self, circuits_to_keep):
         -------
         CircuitList
         """
-        if isinstance(circuits_to_keep, set):
-            new_circuits = list(filter(lambda c: c in circuits_to_keep, self._circuits))
-        else:
-            current_circuits = set(self._circuits)
-            new_circuits = list(filter(lambda c: c in current_circuits, circuits_to_keep))
+        new_circuits = list(filter(lambda c: c in set(circuits_to_keep), self._circuits))
         return CircuitList(new_circuits, self.op_label_aliases)  # don't transfer weights or name
 
     def truncate_to_dataset(self, dataset):
diff --git a/pygsti/protocols/rb.py b/pygsti/protocols/rb.py
index 484e5fb68..9929752fe 100644
--- a/pygsti/protocols/rb.py
+++ b/pygsti/protocols/rb.py
@@ -10,6 +10,7 @@
 # http://www.apache.org/licenses/LICENSE-2.0 or in the LICENSE file in the root pyGSTi directory.
 #***************************************************************************************************
 
+from collections import defaultdict
 import numpy as _np
 
 from pygsti.protocols import protocol as _proto
@@ -111,7 +112,7 @@ class CliffordRBDesign(_vb.BenchmarkingDesign):
     """
 
     @classmethod
-    def from_existing_circuits(cls, circuits_and_idealouts_by_depth, qubit_labels=None,
+    def from_existing_circuits(cls, data_by_depth, qubit_labels=None,
                                randomizeout=False, citerations=20, compilerargs=(), interleaved_circuit=None,
                                descriptor='A Clifford RB experiment', add_default_protocol=False):
         """
@@ -123,10 +124,12 @@ def from_existing_circuits(cls, circuits_and_idealouts_by_depth, qubit_labels=No
 
         Parameters
         ----------
-        circuits_and_idealouts_by_depth : dict
-            A dictionary whose keys are integer depths and whose values are lists of `(circuit, ideal_outcome)` 
-            2-tuples giving each RB circuit and its 
-            ideal (correct) outcome.
+        data_by_depth : dict
+            A dictionary whose keys are integer depths and whose values are lists of
+            `(circuit, ideal_outcome, num_native_gates)` tuples giving each RB circuit, its 
+            ideal (correct) outcome, and (optionally) the number of native gates in the compiled Cliffords.
+            If only a 2-tuple is passed, i.e. number of native gates is not included,
+            the :meth:`average_gates_per_clifford()` function will not work.
 
         qubit_labels : list, optional
             If not None, a list of the qubits that the RB circuits are to be sampled for. This should
@@ -183,22 +186,27 @@ def from_existing_circuits(cls, circuits_and_idealouts_by_depth, qubit_labels=No
         -------
         CliffordRBDesign
         """
-        depths = sorted(list(circuits_and_idealouts_by_depth.keys()))
-        circuit_lists = [[x[0] for x in circuits_and_idealouts_by_depth[d]] for d in depths]
-        ideal_outs = [[x[1] for x in circuits_and_idealouts_by_depth[d]] for d in depths]
-        circuits_per_depth = [len(circuits_and_idealouts_by_depth[d]) for d in depths]
+        depths = sorted(list(data_by_depth.keys()))
+        circuit_lists = [[x[0] for x in data_by_depth[d]] for d in depths]
+        ideal_outs = [[x[1] for x in data_by_depth[d]] for d in depths]
+        try:
+            native_gate_counts = [[x[2] for x in data_by_depth[d]] for d in depths]
+        except IndexError:
+            native_gate_counts = None
+        circuits_per_depth = [len(data_by_depth[d]) for d in depths]
         self = cls.__new__(cls)
         self._init_foundation(depths, circuit_lists, ideal_outs, circuits_per_depth, qubit_labels,
                               randomizeout, citerations, compilerargs, descriptor, add_default_protocol,
-                              interleaved_circuit)
+                              interleaved_circuit, native_gate_counts=native_gate_counts)
         return self
 
     def __init__(self, pspec, clifford_compilations, depths, circuits_per_depth, qubit_labels=None, randomizeout=False,
-                 interleaved_circuit=None, citerations=20, compilerargs=(), descriptor='A Clifford RB experiment',
-                 add_default_protocol=False, seed=None, verbosity=1, num_processes=1):
+                 interleaved_circuit=None, citerations=20, compilerargs=(), exact_compilation_key=None,
+                 descriptor='A Clifford RB experiment', add_default_protocol=False, seed=None, verbosity=1, num_processes=1):
         if qubit_labels is None: qubit_labels = tuple(pspec.qubit_labels)
         circuit_lists = []
         ideal_outs = []
+        native_gate_counts = []
 
         if seed is None:
             self.seed = _np.random.randint(1, 1e6)  # Pick a random seed
@@ -214,26 +222,35 @@ def __init__(self, pspec, clifford_compilations, depths, circuits_per_depth, qub
             args_list = [(pspec, clifford_compilations, l)] * circuits_per_depth
             kwargs_list = [dict(qubit_labels=qubit_labels, randomizeout=randomizeout, citerations=citerations,
                                 compilerargs=compilerargs, interleaved_circuit=interleaved_circuit,
-                                seed=lseed + i) for i in range(circuits_per_depth)]
+                                seed=lseed + i, return_native_gate_counts=True, exact_compilation_key=exact_compilation_key)
+                                for i in range(circuits_per_depth)]
             results = _tools.mptools.starmap_with_kwargs(_rc.create_clifford_rb_circuit, circuits_per_depth,
                                                          num_processes, args_list, kwargs_list)
 
             circuits_at_depth = []
             idealouts_at_depth = []
-            for c, iout in results:
+            native_gate_counts_at_depth = []
+            for c, iout, nng in results:
                 circuits_at_depth.append(c)
                 idealouts_at_depth.append((''.join(map(str, iout)),))
+                native_gate_counts_at_depth.append(nng)
 
             circuit_lists.append(circuits_at_depth)
             ideal_outs.append(idealouts_at_depth)
+            native_gate_counts.append(native_gate_counts_at_depth)
 
         self._init_foundation(depths, circuit_lists, ideal_outs, circuits_per_depth, qubit_labels,
                               randomizeout, citerations, compilerargs, descriptor, add_default_protocol,
-                              interleaved_circuit)
+                              interleaved_circuit, native_gate_counts=native_gate_counts)
 
     def _init_foundation(self, depths, circuit_lists, ideal_outs, circuits_per_depth, qubit_labels,
                          randomizeout, citerations, compilerargs, descriptor, add_default_protocol,
-                         interleaved_circuit):
+                         interleaved_circuit, native_gate_counts=None, exact_compilation_key=None):
+        self.native_gate_count_lists = native_gate_counts
+        if self.native_gate_count_lists is not None:
+            # If we have native gate information, pair this with circuit data so that we serialize/truncate properly
+            self.paired_with_circuit_attrs = ["native_gate_count_lists"]
+
         super().__init__(depths, circuit_lists, ideal_outs, qubit_labels, remove_duplicates=False)
         self.circuits_per_depth = circuits_per_depth
         self.randomizeout = randomizeout
@@ -241,6 +258,7 @@ def _init_foundation(self, depths, circuit_lists, ideal_outs, circuits_per_depth
         self.compilerargs = compilerargs
         self.descriptor = descriptor
         self.interleaved_circuit = interleaved_circuit
+        self.exact_compilation_key = exact_compilation_key
         if add_default_protocol:
             if randomizeout:
                 defaultfit = 'A-fixed'
@@ -248,6 +266,92 @@ def _init_foundation(self, depths, circuit_lists, ideal_outs, circuits_per_depth
                 defaultfit = 'full'
             self.add_default_protocol(RB(name='RB', defaultfit=defaultfit))
 
+    def average_native_gates_per_clifford_for_circuit(self, list_idx, circ_idx):
+        """The average number of native gates per Clifford for a specific circuit
+
+        Parameters
+        ----------
+        list_idx: int
+            The index of the circuit list (for a given depth)
+
+        circ_idx: int
+            The index of the circuit within the circuit list
+        
+        Returns
+        -------
+        avg_gate_counts: dict
+            The average number of native gates, native 2Q gates, and native size
+            per Clifford as values with respective label keys
+        """
+        if self.native_gate_counts_lists is None:
+            raise ValueError("Native gate counts not available, cannot compute average gates per Clifford")
+        
+        num_clifford_gates = self.depths[list_idx] + 1
+        avg_gate_counts = {}
+        for key, native_gate_count in self.native_gate_count_lists[list_idx][circ_idx].items():
+            avg_gate_counts[key.replace('native', 'avg_native_per_clifford')] = native_gate_count / num_clifford_gates
+        
+        return avg_gate_counts
+
+    def average_native_gates_per_clifford_for_circuit_list(self, list_idx):
+        """The average number of gates per Clifford for a circuit list
+
+        This essentially gives the average number of native gates per Clifford
+        for a given depth (indexed by list index, not depth).
+
+        Parameters
+        ----------
+        list_idx: int
+            The index of the circuit list (for a given depth)
+
+        circ_idx: int
+            The index of the circuit within the circuit list
+        
+        Returns
+        -------
+        float
+            The average number of native gates per Clifford
+        """
+        if self.native_gate_count_lists is None:
+            raise ValueError("Native gate counts not available, cannot compute average gates per Clifford")
+        
+        gate_counts = defaultdict(int)
+        for native_gate_counts in self.native_gate_count_lists[list_idx]:
+            for k, v in native_gate_counts.items():
+                gate_counts[k] += v
+        
+        num_clifford_gates = len(self.native_gate_count_lists[list_idx]) * (self.depths[list_idx] + 1)
+        avg_gate_counts = {}
+        for key, total_native_gate_counts in gate_counts.items():
+            avg_gate_counts[key.replace('native', 'avg_native_per_clifford')] = total_native_gate_counts / num_clifford_gates
+        
+        return avg_gate_counts
+
+    def average_native_gates_per_clifford(self):
+        """The average number of native gates per Clifford for all circuits
+
+        Returns
+        -------
+        float
+            The average number of native gates per Clifford
+        """
+        if self.native_gate_count_lists is None:
+            raise ValueError("Number of native gates not available, cannot compute average gates per Clifford")
+        
+        gate_counts = defaultdict(int)
+        num_clifford_gates = 0
+        for list_idx in range(len(self.depths)):
+            for native_gate_counts in self.native_gate_count_lists[list_idx]:
+                for k, v in native_gate_counts.items():
+                    gate_counts[k] += v
+            num_clifford_gates += len(self.native_gate_count_lists[list_idx]) * (self.depths[list_idx] + 1)
+            
+        avg_gate_counts = {}
+        for key, total_native_gate_counts in gate_counts.items():
+            avg_gate_counts[key.replace('native', 'avg_native_per_clifford')] = total_native_gate_counts / num_clifford_gates
+
+        return avg_gate_counts
+
     def map_qubit_labels(self, mapper):
         """
         Creates a new experiment design whose circuits' qubit labels are updated according to a given mapping.
@@ -970,6 +1074,9 @@ def __init__(self, pspec, clifford_compilations, depths, circuits_per_depth, qub
     def _init_foundation(self, depths, circuit_lists, measurements, signs, circuits_per_depth, qubit_labels, layer_sampling,
                          sampler, samplerargs,  addlocal, lsargs, descriptor,
                          add_default_protocol):
+        # Pair these attributes with circuit data so that we serialize/truncate properly
+        self.paired_with_circuit_attrs = ["measurements", "signs"]
+
         super().__init__(depths, circuit_lists, signs, qubit_labels, remove_duplicates=False)
         self.measurements = measurements
         self.signs = signs
@@ -986,10 +1093,6 @@ def _init_foundation(self, depths, circuit_lists, measurements, signs, circuits_
 
         defaultfit = 'A-fixed'
         self.add_default_protocol(RB(name='RB', defaultfit=defaultfit))
-            
-        self.auxfile_types['signs'] = 'json' # Makes sure that signs and measurements are saved seperately
-        self.auxfile_types['measurements'] = 'json'
-
 
 class RandomizedBenchmarking(_vb.SummaryStatistics):
     """
diff --git a/pygsti/protocols/vb.py b/pygsti/protocols/vb.py
index 7c52e81a6..ea20b64a6 100644
--- a/pygsti/protocols/vb.py
+++ b/pygsti/protocols/vb.py
@@ -11,11 +11,12 @@
 #***************************************************************************************************
 
 import numpy as _np
+import copy as _copy
 
-from pygsti.protocols import protocol as _proto
-from pygsti.models.oplessmodel import SuccessFailModel as _SuccessFailModel
 from pygsti import tools as _tools
 from pygsti.algorithms import randomcircuit as _rc
+from pygsti.protocols import protocol as _proto
+from pygsti.models.oplessmodel import SuccessFailModel as _SuccessFailModel
 
 
 class ByDepthDesign(_proto.CircuitListsDesign):
@@ -67,6 +68,25 @@ def map_qubit_labels(self, mapper):
         mapped_qubit_labels = self._mapped_qubit_labels(mapper)
         return ByDepthDesign(self.depths, mapped_circuit_lists, mapped_qubit_labels, remove_duplicates=False)
 
+    def truncate_to_lists(self, list_indices_to_keep):
+        """
+        Truncates this experiment design by only keeping a subset of its circuit lists.
+
+        Parameters
+        ----------
+        list_indices_to_keep : iterable
+            A list of the (integer) list indices to keep.
+
+        Returns
+        -------
+        ByDepthDesign
+            The truncated experiment design.
+        """
+        ret = _copy.deepcopy(self) # Works for derived classes too
+        ret.depths = [self.depths[i] for i in list_indices_to_keep]
+        ret.circuit_lists = [self.circuit_lists[i] for i in list_indices_to_keep]
+        return ret
+
 
 class BenchmarkingDesign(ByDepthDesign):
     """
@@ -98,12 +118,27 @@ class BenchmarkingDesign(ByDepthDesign):
         Whether to remove duplicates when automatically creating
         all the circuits that need data.
     """
+    
+    paired_with_circuit_attrs = None
+    """List of attributes which are paired up with circuit lists
+
+    These will be saved as external files during serialization,
+    and are truncated when circuit lists are truncated.
+    """
 
     def __init__(self, depths, circuit_lists, ideal_outs, qubit_labels=None, remove_duplicates=False):
         assert(len(depths) == len(ideal_outs))
         super().__init__(depths, circuit_lists, qubit_labels, remove_duplicates)
+        
         self.idealout_lists = ideal_outs
-        self.auxfile_types['idealout_lists'] = 'json'
+
+        if self.paired_with_circuit_attrs is None:
+            self.paired_with_circuit_attrs = ['idealout_lists']
+        else:
+            self.paired_with_circuit_attrs.insert(0, 'idealout_lists')
+        
+        for paired_attr in self.paired_with_circuit_attrs:
+            self.auxfile_types[paired_attr] = 'json'
 
     def _mapped_circuits_and_idealouts_by_depth(self, mapper):
         """ Used in derived classes """
@@ -133,6 +168,76 @@ def map_qubit_labels(self, mapper):
         mapped_qubit_labels = self._mapped_qubit_labels(mapper)
         return BenchmarkingDesign(self.depths, mapped_circuit_lists, list(self.idealout_lists),
                                   mapped_qubit_labels, remove_duplicates=False)
+    
+    def truncate_to_lists(self, list_indices_to_keep):
+        """
+        Truncates this experiment design by only keeping a subset of its circuit lists.
+
+        Parameters
+        ----------
+        list_indices_to_keep : iterable
+            A list of the (integer) list indices to keep.
+
+        Returns
+        -------
+        BenchmarkingDesign
+            The truncated experiment design.
+        """
+        ret = _copy.deepcopy(self) # Works for derived classes too
+        ret.depths = [self.depths[i] for i in list_indices_to_keep]
+        ret.circuit_lists = [self.circuit_lists[i] for i in list_indices_to_keep]
+        for paired_attr in self.paired_with_circuit_attrs:
+            val = getattr(self, paired_attr)
+            new_val = [val[i] for i in list_indices_to_keep]
+            setattr(ret, paired_attr, new_val)
+        return ret
+
+    def _truncate_to_circuits_inplace(self, circuits_to_keep):
+        truncated_circuit_lists = []
+        paired_attr_lists_list = [getattr(self, paired_attr) for paired_attr in self.paired_with_circuit_attrs]
+        truncated_paired_attr_lists_list = [[] for _ in range(len(self.paired_with_circuit_attrs))]
+        for list_idx, circuits in enumerate(self.circuit_lists):
+            paired_attrs = [pal[list_idx] for pal in paired_attr_lists_list]
+            # Do the same filtering as CircuitList.truncate, but drag along any paired attributes
+            new_data = list(zip(*filter(lambda ci: ci[0] in set(circuits_to_keep), zip(circuits, *paired_attrs))))
+            if len(new_data):
+                truncated_circuit_lists.append(new_data[0])
+                for i, attr_data in enumerate(new_data[1:]):
+                    truncated_paired_attr_lists_list[i].append(attr_data)
+            else:
+                # If we have truncated all circuits, append empty lists
+                truncated_circuit_lists.append([])
+                truncated_paired_attr_lists_list.append([[] for _ in range(len(self.paired_with_circuit_attrs))])
+
+        self.circuit_lists = truncated_circuit_lists
+        for paired_attr, paired_attr_lists in zip(self.paired_with_circuit_attrs, truncated_paired_attr_lists_list):
+            setattr(self, paired_attr, paired_attr_lists)
+        super()._truncate_to_circuits_inplace(circuits_to_keep)
+
+    def _truncate_to_design_inplace(self, other_design):
+        truncated_circuit_lists = []
+        paired_attr_lists_list = [getattr(self, paired_attr) for paired_attr in self.paired_with_circuit_attrs]
+        truncated_paired_attr_lists_list = [[] for _ in range(len(self.paired_with_circuit_attrs))]
+        for list_idx, circuits in enumerate(self.circuit_lists):
+            paired_attrs = [pal[list_idx] for pal in paired_attr_lists_list]
+            # Do the same filtering as CircuitList.truncate, but drag along any paired attributes
+            new_data = list(zip(*filter(lambda ci: ci[0] in set(other_design.circuit_lists[list_idx]), zip(circuits, *paired_attrs))))
+            if len(new_data):
+                truncated_circuit_lists.append(new_data[0])
+                for i, attr_data in enumerate(new_data[1:]):
+                    truncated_paired_attr_lists_list[i].append(attr_data)
+            else:
+                # If we have truncated all circuits, append empty lists
+                truncated_circuit_lists.append([])
+                truncated_paired_attr_lists_list.append([[] for _ in range(len(self.paired_with_circuit_attrs))])
+
+        self.circuit_lists = truncated_circuit_lists
+        for paired_attr, paired_attr_lists in zip(self.paired_with_circuit_attrs, truncated_paired_attr_lists_list):
+            setattr(self, paired_attr, paired_attr_lists)
+        super()._truncate_to_design_inplace(other_design)
+
+    def _truncate_to_available_data_inplace(self, dataset):
+        self._truncate_to_circuits_inplace(set(dataset.keys()))
 
 
 class PeriodicMirrorCircuitDesign(BenchmarkingDesign):
diff --git a/test/unit/protocols/test_protocols.py b/test/unit/protocols/test_protocols.py
index 1385707b6..1a5d3f4bc 100644
--- a/test/unit/protocols/test_protocols.py
+++ b/test/unit/protocols/test_protocols.py
@@ -11,6 +11,51 @@ class ExperimentDesignTester(BaseCase):
     def setUpClass(cls):
         cls.gst_design = std.create_gst_experiment_design(4)
 
+        #Create a bunch of experiment designs:
+        from pygsti.protocols import ExperimentDesign, CircuitListsDesign, CombinedExperimentDesign, \
+            SimultaneousExperimentDesign, FreeformDesign, StandardGSTDesign, GateSetTomographyDesign, \
+            CliffordRBDesign, DirectRBDesign, MirrorRBDesign
+        from pygsti.processors import CliffordCompilationRules as CCR
+
+        circuits_on0 = pygsti.circuits.to_circuits(["{}@(0)", "Gxpi2:0", "Gypi2:0"], line_labels=(0,))
+        circuits_on0b = pygsti.circuits.to_circuits(["Gxpi2:0^2", "Gypi2:0^2"], line_labels=(0,))
+        circuits_on1 = pygsti.circuits.to_circuits(["Gxpi2:1^2", "Gypi2:1^2"], line_labels=(1,))
+        circuits_on01 = pygsti.circuits.to_circuits(["Gcnot:0:1", "Gxpi2:0Gypi2:1^2Gcnot:0:1Gxpi:0"],
+                                                    line_labels=(0,1))
+
+        #For GST edesigns
+        mdl = std.target_model()
+        gst_pspec = mdl.create_processor_spec()
+
+        #For RB edesigns
+        pspec = pygsti.processors.QubitProcessorSpec(2, ["Gxpi2", "Gypi2","Gxx"],
+                                                     geometry='line', qubit_labels=(0,1))
+        compilations = {"absolute": CCR.create_standard(pspec, "absolute", ("paulis", "1Qcliffords"), verbosity=0),
+                        "paulieq": CCR.create_standard(pspec, "paulieq", ("1Qcliffords", "allcnots"), verbosity=0),
+                        }
+
+        pspec1Q = pygsti.processors.QubitProcessorSpec(1, ["Gxpi2", "Gypi2","Gxmpi2", "Gympi2"],
+                                                       geometry='line', qubit_labels=(0,))
+        compilations1Q = {"absolute": CCR.create_standard(pspec1Q, "absolute", ("paulis", "1Qcliffords"), verbosity=0),
+                          "paulieq": CCR.create_standard(pspec1Q, "paulieq", ("1Qcliffords", "allcnots"), verbosity=0),
+                          }
+
+        edesigns = []
+        edesigns.append(ExperimentDesign(circuits_on0))
+        edesigns.append(CircuitListsDesign([circuits_on0, circuits_on0b]))
+        edesigns.append(CombinedExperimentDesign({'one': ExperimentDesign(circuits_on0),
+                                                  'two': ExperimentDesign(circuits_on1),
+                                                  'three': ExperimentDesign(circuits_on01)}, qubit_labels=(0,1)))
+        edesigns.append(SimultaneousExperimentDesign([ExperimentDesign(circuits_on0), ExperimentDesign(circuits_on1)]))
+        edesigns.append(FreeformDesign(circuits_on01))
+        edesigns.append(std.create_gst_experiment_design(2))
+        edesigns.append(GateSetTomographyDesign(gst_pspec, [circuits_on0, circuits_on0b]))
+        edesigns.append(CliffordRBDesign(pspec, compilations, depths=[0,2,5], circuits_per_depth=4))
+        edesigns.append(DirectRBDesign(pspec, compilations, depths=[0,2,5], circuits_per_depth=4))
+        edesigns.append(MirrorRBDesign(pspec1Q, depths=[0,2,4], circuits_per_depth=4,
+                                       clifford_compilations=compilations1Q))
+        cls.edesigns = edesigns
+
     def test_promotion(self):
         circuits = pygsti.circuits.to_circuits(["{}@(0)", "Gxpi2:0", "Gypi2:0"])
         edesign1 = pygsti.protocols.ExperimentDesign(circuits)
@@ -89,52 +134,7 @@ def test_create_edesign_fromdir_subdirs(self, root_path):
         self.assertTrue(all([a == b for a,b in zip(edesign3['subdir2'].all_circuits_needing_data, self.gst_design.circuit_lists[1])]))
 
     def test_map_edesign_sslbls(self):
-        #Create a bunch of experiment designs:
-        from pygsti.protocols import ExperimentDesign, CircuitListsDesign, CombinedExperimentDesign, \
-            SimultaneousExperimentDesign, FreeformDesign, StandardGSTDesign, GateSetTomographyDesign, \
-            CliffordRBDesign, DirectRBDesign, MirrorRBDesign
-        from pygsti.processors import CliffordCompilationRules as CCR
-
-        circuits_on0 = pygsti.circuits.to_circuits(["{}@(0)", "Gxpi2:0", "Gypi2:0"], line_labels=(0,))
-        circuits_on0b = pygsti.circuits.to_circuits(["Gxpi2:0^2", "Gypi2:0^2"], line_labels=(0,))
-        circuits_on1 = pygsti.circuits.to_circuits(["Gxpi2:1^2", "Gypi2:1^2"], line_labels=(1,))
-        circuits_on01 = pygsti.circuits.to_circuits(["Gcnot:0:1", "Gxpi2:0Gypi2:1^2Gcnot:0:1Gxpi:0"],
-                                                    line_labels=(0,1))
-
-        #For GST edesigns
-        mdl = std.target_model()
-        gst_pspec = mdl.create_processor_spec()
-
-        #For RB edesigns
-        pspec = pygsti.processors.QubitProcessorSpec(2, ["Gxpi2", "Gypi2","Gxx"],
-                                                     geometry='line', qubit_labels=(0,1))
-        compilations = {"absolute": CCR.create_standard(pspec, "absolute", ("paulis", "1Qcliffords"), verbosity=0),
-                        "paulieq": CCR.create_standard(pspec, "paulieq", ("1Qcliffords", "allcnots"), verbosity=0),
-                        }
-
-        pspec1Q = pygsti.processors.QubitProcessorSpec(1, ["Gxpi2", "Gypi2","Gxmpi2", "Gympi2"],
-                                                       geometry='line', qubit_labels=(0,))
-        compilations1Q = {"absolute": CCR.create_standard(pspec1Q, "absolute", ("paulis", "1Qcliffords"), verbosity=0),
-                          "paulieq": CCR.create_standard(pspec1Q, "paulieq", ("1Qcliffords", "allcnots"), verbosity=0),
-                          }
-
-
-        edesigns = []
-        edesigns.append(ExperimentDesign(circuits_on0))
-        edesigns.append(CircuitListsDesign([circuits_on0, circuits_on0b]))
-        edesigns.append(CombinedExperimentDesign({'one': ExperimentDesign(circuits_on0),
-                                                  'two': ExperimentDesign(circuits_on1),
-                                                  'three': ExperimentDesign(circuits_on01)}, qubit_labels=(0,1)))
-        edesigns.append(SimultaneousExperimentDesign([ExperimentDesign(circuits_on0), ExperimentDesign(circuits_on1)]))
-        edesigns.append(FreeformDesign(circuits_on01))
-        edesigns.append(std.create_gst_experiment_design(2))
-        edesigns.append(GateSetTomographyDesign(gst_pspec, [circuits_on0, circuits_on0b]))
-        edesigns.append(CliffordRBDesign(pspec, compilations, depths=[0,2,5], circuits_per_depth=4))
-        edesigns.append(DirectRBDesign(pspec, compilations, depths=[0,2,5], circuits_per_depth=4))
-        edesigns.append(MirrorRBDesign(pspec1Q, depths=[0,2,4], circuits_per_depth=4,
-                                       clifford_compilations=compilations1Q))
-
-        for edesign in edesigns:
+        for edesign in self.edesigns:
             print("Testing edesign of type: ", str(type(edesign)))
             orig_qubits = edesign.qubit_labels
             for c in edesign.all_circuits_needing_data:
@@ -150,3 +150,28 @@ def test_map_edesign_sslbls(self):
             self.assertEqual(mapped_edesign.qubit_labels, mapped_qubits)
             for c in mapped_edesign.all_circuits_needing_data:
                 self.assertTrue(set(c.line_labels).issubset(mapped_qubits))
+    
+    def test_truncation(self):
+        from pygsti.protocols import BenchmarkingDesign
+        
+        for edesign in self.edesigns:
+            print("Testing edesign of type: ", str(type(edesign)))
+            
+            truncated_circuits = edesign.all_circuits_needing_data[:2]
+            truncated_edesign = edesign.truncate_to_circuits(truncated_circuits)
+            self.assertTrue(set(truncated_circuits) == set(truncated_edesign.all_circuits_needing_data))
+
+            if isinstance(edesign, BenchmarkingDesign):
+                # Check that the paired attributes were also truncated properly
+                # These will be lists of lists
+                for attr in edesign.paired_with_circuit_attrs:
+                    attr_lists = getattr(edesign, attr)
+                    truncated_attr_lists = getattr(truncated_edesign, attr)
+                    for a_list, ta_list, c_list, tc_list in zip(attr_lists, truncated_attr_lists,
+                                                                edesign.circuit_lists, truncated_edesign.circuit_lists):
+                        self.assertTrue(len(ta_list) == len(tc_list))
+
+                        # Ensure that the paired attribute data is correct
+                        for ta_data, tc_circ in zip(ta_list, tc_list):
+                            untruncated_idx = c_list.index(tc_circ)
+                            self.assertTrue(a_list[untruncated_idx] == ta_data)
diff --git a/test/unit/protocols/test_rb.py b/test/unit/protocols/test_rb.py
index 7cf601b7c..19259e5a2 100644
--- a/test/unit/protocols/test_rb.py
+++ b/test/unit/protocols/test_rb.py
@@ -1,5 +1,7 @@
 from ..util import BaseCase
 
+import numpy as _np
+
 import pygsti
 from pygsti.protocols import rb as _rb
 from pygsti.processors import CliffordCompilationRules as CCR
@@ -11,7 +13,7 @@ def setUp(self):
         self.num_qubits = 2
         self.qubit_labels = ['Q'+str(i) for i in range(self.num_qubits)]
         
-        gate_names = ['Gxpi2', 'Gxmpi2', 'Gypi2', 'Gympi2', 'Gcphase']
+        gate_names = ['Gi', 'Gxpi2', 'Gxmpi2', 'Gypi2', 'Gympi2', 'Gcphase']
         availability = {'Gcphase':[('Q'+str(i),'Q'+str((i+1) % self.num_qubits)) for i in range(self.num_qubits)]}
 
         self.pspec = pygsti.processors.QubitProcessorSpec(self.num_qubits, gate_names, availability=availability,
@@ -21,7 +23,16 @@ def setUp(self):
             'paulieq': CCR.create_standard(self.pspec, 'paulieq', ('1Qcliffords', 'allcnots'), verbosity=0)
         }
 
+        gate_names_1Q = gate_names[:-1]
+        self.qubit_labels1Q = ['Q0']
+        self.pspec1Q = pygsti.processors.QubitProcessorSpec(1, gate_names_1Q, qubit_labels=self.qubit_labels1Q)
+        self.compilations1Q = {
+            'absolute': CCR.create_standard(self.pspec1Q, 'absolute', ('paulis', '1Qcliffords'), verbosity=0),
+            'paulieq': CCR.create_standard(self.pspec1Q, 'paulieq', ('1Qcliffords', 'allcnots'), verbosity=0)
+        }
+
         # TODO: Test a lot of these, currently just the default from the tutorial
+        # Probably as pytest mark parameterize for randomizeout, compilerargs?
         self.depths = [0, 2]#, 4, 8]
         self.circuits_per_depth = 5
         self.qubits = ['Q0', 'Q1']
@@ -61,7 +72,39 @@ def test_design_construction(self):
 
         [[self.assertAlmostEqual(c.simulate(tmodel)[bs],1.) for c, bs in zip(cl, bsl)] for cl, bsl in zip(mp_design.circuit_lists, mp_design.idealout_lists)]
 
+    def test_deterministic_compilation(self):        
+        # TODO: Figure out good test for this. Full circuit is a synthetic idle, we need to somehow check the non-inverted
+        # Clifford is the same as the random case?
+        abs_design = _rb.CliffordRBDesign(
+            self.pspec1Q, self.compilations1Q, self.depths, self.circuits_per_depth, qubit_labels=self.qubit_labels1Q,
+            randomizeout=self.randomizeout, interleaved_circuit=self.interleaved_circuit,
+            citerations=self.citerations, compilerargs=self.compiler_args, seed=self.seed,
+            verbosity=self.verbosity, exact_compilation_key='absolute')
+        
+        peq_design = _rb.CliffordRBDesign(
+            self.pspec1Q, self.compilations1Q, self.depths, self.circuits_per_depth, qubit_labels=self.qubit_labels1Q,
+            randomizeout=self.randomizeout, interleaved_circuit=self.interleaved_circuit,
+            citerations=self.citerations, compilerargs=self.compiler_args, seed=self.seed,
+            verbosity=self.verbosity, exact_compilation_key='paulieq')
+
+        # Testing a non-standard (but unrealistic) compilation
+        rule_dict = {f'C{i}': (_np.eye(2), pygsti.circuits.Circuit([], (0,))) for i in range(24)}
+        compilations = self.compilations1Q.copy()
+        compilations["idle"] = pygsti.processors.CompilationRules(rule_dict)
+        idle_design = _rb.CliffordRBDesign(
+            self.pspec1Q, compilations, self.depths, self.circuits_per_depth, qubit_labels=self.qubit_labels1Q,
+            randomizeout=False, interleaved_circuit=self.interleaved_circuit,
+            citerations=self.citerations, compilerargs=self.compiler_args, seed=self.seed,
+            verbosity=self.verbosity, exact_compilation_key='idle')
+
+        # All circuits should be the empty circuit (since we've turned off randomizeout)
+        for clist in idle_design.circuit_lists:
+            self.assertTrue(set(clist) == set([pygsti.circuits.Circuit([], self.qubit_labels1Q)]))
 
+        # Also a handy place to test native gate counts since it should be 0
+        avg_gate_counts = idle_design.average_native_gates_per_clifford()
+        for v in avg_gate_counts.values():
+            self.assertTrue(v == 0)
 
 class TestDirectRBDesign(BaseCase):