Perform exhaustive symmetrization

    :param qc: A QuantumComputer which can run quantum programs
    :param qubits: qubits on which the symmetrization program runs
    :param shots: number of shots in the symmetrized program
    :prog: program to symmetrize
    :return: - the equivalent of a `run` output, but with exhaustive symmetrization
             - dict keyed by qubit, valued by index of the numpy array containing
                    bitstring results

def _exhaustive_symmetrization(qc: QuantumComputer, qubits: List[int],
                               shots: int, prog: Program) -> (np.ndarray, Dict):
    """
    Perform exhaustive symmetrization

    :param qc: A QuantumComputer which can run quantum programs
    :param qubits: qubits on which the symmetrization program runs
    :param shots: number of shots in the symmetrized program
    :prog: program to symmetrize
    :return: - the equivalent of a `run` output, but with exhaustive symmetrization
             - dict keyed by qubit, valued by index of the numpy array containing
                    bitstring results
    """
    # Symmetrize -- flip qubits pre-measurement
    n_shots_symm = int(round(np.ceil(shots / 2**len(qubits))))
    if n_shots_symm * 2**len(qubits) > shots:
        warnings.warn(f"Symmetrization increasing number of shots from {shots} to {round(n_shots_symm * 2**len(qubits))}")
    list_bitstrings_symm = []
    for ops_bool in itertools.product([0, 1], repeat=len(qubits)):
        total_prog_symm = prog.copy()
        prog_symm = _ops_bool_to_prog(ops_bool, qubits)
        total_prog_symm += prog_symm
        # Run the experiment
        dict_qub_idx = {}
        ro = total_prog_symm.declare('ro', 'BIT', len(qubits))
        for i, q in enumerate(qubits):
            total_prog_symm += MEASURE(q, ro[i])
            # Keep track of qubit-classical register mapping via dict
            dict_qub_idx[q] = i
        total_prog_symm.wrap_in_numshots_loop(n_shots_symm)
        total_prog_symm_native = qc.compiler.quil_to_native_quil(total_prog_symm)
        total_prog_symm_bin = qc.compiler.native_quil_to_executable(total_prog_symm_native)
        bitstrings_symm = qc.run(total_prog_symm_bin)
        # Flip the results post-measurement
        bitstrings_symm = bitstrings_symm ^ ops_bool
        # Gather together the symmetrized results into list
        list_bitstrings_symm.append(bitstrings_symm)

    # Gather together all the symmetrized results
    bitstrings = reduce(lambda x, y: np.vstack((x, y)), list_bitstrings_symm)
    return bitstrings, dict_qub_idx