diff --git a/symmer/operators/base.py b/symmer/operators/base.py
index 8096cfe8..80f87760 100644
--- a/symmer/operators/base.py
+++ b/symmer/operators/base.py
@@ -1104,44 +1104,7 @@ def is_noncontextual(self) -> bool:
         if self.n_terms < 4:
             # all operators with 3 or less P are noncontextual
             return True
-
-        to_reduce = np.vstack([np.hstack([self.Z_block, self.X_block]), np.eye(2 * self.n_qubits, dtype=bool)])
-        cref_matrix = _cref_binary(to_reduce)
-        Z2_symp = cref_matrix[self.n_terms:, np.all(~cref_matrix[:self.n_terms], axis=0)].T
-        Z2_terms = PauliwordOp(Z2_symp, np.ones(Z2_symp.shape[0]))
-
-        if Z2_terms.n_terms < 1:
-            remaining = self.cleanup()
-            if remaining.n_terms > 2 * remaining.n_qubits + 1:
-                # no symmetry component, therefore operator must be made up of pairwise
-                # anticommuting pauli operators to be noncontextual. This can have a max size of 2n+1, if larger cannot
-                # be noncontextual
-                return False
-            else:
-                # for remaining to be noncontextual must be disjoint union of cliques
-                # we can test for this below
-                adj_matrix_view = np.ascontiguousarray(remaining.adjacency_matrix).view(
-                    np.dtype((np.void, remaining.adjacency_matrix.dtype.itemsize * remaining.adjacency_matrix.shape[1]))
-                )
-                re_order_indices = np.argsort(adj_matrix_view.ravel())
-                # sort the adj matrix and vector of coefficients accordingly
-                sorted_terms = remaining.adjacency_matrix[re_order_indices]
-                # unique terms are those with non-zero entries in the adjacent row difference array
-                diff_adjacent = np.diff(sorted_terms, axis=0)
-                mask_unique_terms = np.append(True, np.any(diff_adjacent, axis=1))
-                clique_mask = sorted_terms[mask_unique_terms]
-
-                # check for overlap (array of ones == no overlap)
-                return np.all(np.sum(clique_mask, axis=0) == 1)
-        else:
-            gens = self.generators
-            if check_adjmat_noncontextual(gens.adjacency_matrix):
-                return True
-            from symmer.utils import get_generators_including_xz_products
-            gens_xyz = get_generators_including_xz_products(self)
-            return check_adjmat_noncontextual(gens_xyz.adjacency_matrix)
-
-
+        return check_adjmat_noncontextual(self.adjacency_matrix)
 
     def _rotate_by_single_Pword(self,
             Pword: "PauliwordOp", 
diff --git a/symmer/utils.py b/symmer/utils.py
index 98242b63..92eb4209 100644
--- a/symmer/utils.py
+++ b/symmer/utils.py
@@ -416,39 +416,4 @@ def Get_AC_root(power: float, operator: AntiCommutingOp) -> PauliwordOp:
 
     AC_root = (rot.dagger * Ps_root * rot).multiply_by_constant(gamma_l ** power)
 
-    return AC_root
-
-
-def get_generators_including_xz_products(operator: PauliwordOp) -> PauliwordOp:
-    """
-    Given an input operator perform similar process to finding generators (but do not use X,Z terms to generate Y terms)
-    i.e. the set can be dependent.
-
-    Args:
-        operator:
-
-    Returns:
-        PauliwordOp of generators build of X,Y,Z terms (note can have dependent terms!)
-
-    """
-    X_only = np.logical_and(operator.X_block, ~operator.Z_block)
-    Z_only = np.logical_and(~operator.X_block, operator.Z_block)
-    Y_only = np.logical_and(operator.Z_block, operator.X_block)
-
-    XZY_block = np.hstack((np.hstack((X_only, Z_only)), Y_only))
-
-    row_red_XYZ = _rref_binary(XZY_block)
-    non_zero_rows = row_red_XYZ[np.sum(row_red_XYZ, axis=1).astype(bool)]
-
-    X_new = non_zero_rows[:, :operator.X_block.shape[1]]
-    Z_new = non_zero_rows[:, operator.X_block.shape[1]:2 * operator.Z_block.shape[1]]
-    Y_new = non_zero_rows[:, 2 * operator.Z_block.shape[1]:]
-
-    sym_final = np.hstack((np.logical_or(X_new, Y_new),
-                           np.logical_or(Z_new, Y_new)))
-
-    # remove duplicates due to same X_new and Z_new positions
-    xyz_gens = PauliwordOp(sym_final, np.ones(sym_final.shape[0])).cleanup()
-    xyz_gens.coeff_vec = np.ones_like(xyz_gens.coeff_vec)
-
-    return xyz_gens
\ No newline at end of file
+    return AC_root
\ No newline at end of file