diff --git a/symmer/projection/contextual_subspace.py b/symmer/projection/contextual_subspace.py
index 2304ddd..0391a9f 100644
--- a/symmer/projection/contextual_subspace.py
+++ b/symmer/projection/contextual_subspace.py
@@ -29,7 +29,7 @@ def manual_stabilizers(self, S: Union[List[str], StabilizerOp]) -> None:
         """
         """
         if isinstance(S, list):
-            S = StabilizerOp.from_list(S_list)
+            S = StabilizerOp.from_list(S)
         self.n_qubits_in_subspace = self.operator.n_qubits - S.n_terms
         self.stabilizers = S
 
diff --git a/symmer/projection/utils.py b/symmer/projection/utils.py
index 62a5adf..b96405d 100644
--- a/symmer/projection/utils.py
+++ b/symmer/projection/utils.py
@@ -72,7 +72,7 @@ def build_basis_weighting_operator(self):
             np.abs(self.weighting_operator.coeff_vec)
         ).cleanup()
         #X_op = self.weighting_operator
-        self.basis_weighting = X_op.sort(key='magnitude')
+        self.basis_weighting = X_op.sort(by='magnitude')
         self.qubit_positions = np.arange(self.weighting_operator.n_qubits)
         self.term_region = [0,self.basis_weighting.n_terms]
         
diff --git a/symmer/symplectic/base.py b/symmer/symplectic/base.py
index f5c8742..2146af9 100644
--- a/symmer/symplectic/base.py
+++ b/symmer/symplectic/base.py
@@ -293,27 +293,27 @@ def copy(self) -> "PauliwordOp":
         return deepcopy(self)
 
     def sort(self, 
-            by: str = 'decreasing',
-            key:str = 'magnitude'
+            by:  str = 'magnitude',
+            key: str = 'decreasing'
         ) -> "PauliwordOp":
         """
-        Sort the terms by some key, either magnitude, weight X, Y or Z
+        Sort the terms either by magnitude, weight X, Y or Z
         """
-        if key=='magnitude':
+        if by=='magnitude':
             sort_order = np.argsort(-abs(self.coeff_vec))
-        elif key=='weight':
-            sort_order = np.argsort(-np.sum(self.symp_matrix, axis=1))
-        elif key=='Z':
-            sort_order = np.argsort(np.sum((self.n_qubits+1)*self.X_block + self.Z_block, axis=1))
-        elif key=='X':
-            sort_order = np.argsort(np.sum(self.X_block + (self.n_qubits+1)*self.Z_block, axis=1))
-        elif key=='Y':
-            sort_order = np.argsort(np.sum(abs(self.X_block - self.Z_block), axis=1))
+        elif by=='weight':
+            sort_order = np.argsort(-np.sum(self.symp_matrix.astype(int), axis=1))
+        elif by=='Z':
+            sort_order = np.argsort(np.sum((self.n_qubits+1)*self.X_block.astype(int) + self.Z_block.astype(int), axis=1))
+        elif by=='X':
+            sort_order = np.argsort(np.sum(self.X_block.astype(int) + (self.n_qubits+1)*self.Z_block.astype(int), axis=1))
+        elif by=='Y':
+            sort_order = np.argsort(np.sum(abs(self.X_block.astype(int) - self.Z_block.astype(int)), axis=1))
         else:
-            raise ValueError('Only permitted sort key values are magnitude, weight, X, Y or Z')
-        if by=='increasing':
+            raise ValueError('Only permitted sort by values are magnitude, weight, X, Y or Z')
+        if key=='increasing':
             sort_order = sort_order[::-1]
-        elif by!='decreasing':
+        elif key!='decreasing':
             raise ValueError('Only permitted sort by values are increasing or decreasing')
         return PauliwordOp(self.symp_matrix[sort_order], self.coeff_vec[sort_order])
 
@@ -572,6 +572,10 @@ def __getitem__(self,
             if stop is None:
                 stop=self.n_terms
             mask = np.arange(start, stop, key.step)
+        elif isinstance(key, (list, np.ndarray)):
+            mask = np.asarray(key)
+        else:
+            raise ValueError('Unrecognised input, must be an integer, slice, list or np.array')
         
         symp_items = self.symp_matrix[mask]
         coeff_items = self.coeff_vec[mask]
@@ -809,7 +813,7 @@ def clique_cover(self,
             if colouring_interchange is not False:
                 warnings.warn(f'{strategy} is not a graph colouring method, so colouring_interchange flag is ignored')
 
-            sorted_op_list = list(self.sort(by='decreasing', key='magnitude'))
+            sorted_op_list = list(self.sort(by='magnitude', key='decreasing'))
 
             check_dic = {
                 'C': lambda x, y: np.all(x.commutes_termwise(y)),
diff --git a/symmer/symplectic/noncontextual_op.py b/symmer/symplectic/noncontextual_op.py
index 2542316..ec08c7b 100644
--- a/symmer/symplectic/noncontextual_op.py
+++ b/symmer/symplectic/noncontextual_op.py
@@ -68,7 +68,7 @@ def _dfs_noncontextual_op(cls, H: PauliwordOp, runtime=10, strategy='magnitude')
         cutoff time ensures if the number of possibilities is large the function will STOP and not take too long
 
         """
-        operator = H.sort(key='magnitude')
+        operator = H.sort(by='magnitude')
         noncontextual_ops = []
 
         n=0
@@ -102,7 +102,7 @@ def _diag_first_noncontextual_op(cls, H: PauliwordOp):
         """
         noncontextual_operator = cls._diag_noncontextual_op(H)
         # order the remaining terms by coefficient magnitude
-        off_diag_terms = (H - noncontextual_operator).sort(key='magnitude')
+        off_diag_terms = (H - noncontextual_operator).sort(by='magnitude')
         # append terms that do not make the noncontextual_operator contextual!
         for term in off_diag_terms:
             if (noncontextual_operator+term).is_noncontextual:
@@ -117,7 +117,7 @@ def _single_sweep_noncontextual_operator(cls, H, strategy='magnitude'):
         noncontextual_operator = PauliwordOp.empty(H.n_qubits)
         
         if strategy=='magnitude':
-            operator = H.sort(key='magnitude')
+            operator = H.sort(by='magnitude')
         elif strategy=='random':
             order = np.arange(H.n_terms)
             np.random.shuffle(order)