How to use the pyquil.paulis.PauliSum function in pyquil

if store_basis:
            self.states = [
                np.binary_repr(i, width=len(self.qubits))
                for i in range(self.nstates)
            ]

        # Check argument types
        if not isinstance(self.cost_ham, (list, tuple)):
            raise TypeError("cost_ham must be a list of PauliSum objects.")
        if not all([isinstance(x, PauliSum) for x in self.cost_ham]):
            raise TypeError("cost_ham must be a list of PauliSum objects")

        if not isinstance(self.ref_ham, (list, tuple)):
            raise TypeError("ref_ham must be a list of PauliSum objects")
        if not all([isinstance(x, PauliSum) for x in self.ref_ham]):
            raise TypeError("ref_ham must be a list of PauliSum objects")

        if not isinstance(self.ref_state_prep, Program):
            raise TypeError("Please provide a pyQuil Program object "
                            "to generate initial state.")

:param vqe_option: (Optional. Default=None). VQE optional arguments.  If None set to
        ``vqe_option = {'disp': print_fun, 'return_all': True, 'samples': samples}``

    """
    if not isinstance(graph, nx.Graph) and isinstance(graph, list):
        maxcut_graph = nx.Graph()
        for edge in graph:
            maxcut_graph.add_edge(*edge)
        graph = maxcut_graph.copy()

    cost_operators = []
    driver_operators = []
    for i, j in graph.edges():
        cost_operators.append(PauliTerm("Z", i, 0.5)*PauliTerm("Z", j) + PauliTerm("I", 0, -0.5))
    for i in graph.nodes():
        driver_operators.append(PauliSum([PauliTerm("X", i, -1.0)]))

    if connection is None:
        connection = get_qc(f"{len(graph.nodes)}q-qvm")

    if minimizer_kwargs is None:
        minimizer_kwargs = {'method': 'Nelder-Mead',
                            'options': {'ftol': 1.0e-2, 'xtol': 1.0e-2,
                                        'disp': False}}
    if vqe_option is None:
        vqe_option = {'disp': print, 'return_all': True,
                      'samples': samples}

    qaoa_inst = QAOA(connection, list(graph.nodes()), steps=steps, cost_ham=cost_operators,
                     ref_ham=driver_operators, store_basis=True,
                     rand_seed=rand_seed,
                     init_betas=initial_beta,

ref_prog.inst(H(i))
            self.ref_state_prep = ref_prog

        if not isinstance(cost_ham, (list, tuple)):
            raise TypeError("""cost_hamiltonian must be a list of PauliSum
                               objects""")
        if not all([isinstance(x, PauliSum) for x in cost_ham]):
            raise TypeError("""cost_hamiltonian must be a list of PauliSum
                                   objects""")
        else:
            self.cost_ham = cost_ham

        if not isinstance(ref_hamiltonian, (list, tuple)):
            raise TypeError("""cost_hamiltonian must be a list of PauliSum
                               objects""")
        if not all([isinstance(x, PauliSum) for x in ref_hamiltonian]):
            raise TypeError("""cost_hamiltonian must be a list of PauliSum
                                   objects""")
        else:
            self.ref_ham = ref_hamiltonian

        if minimizer is None:
            self.minimizer = optimize.minimize
        else:
            self.minimizer = minimizer
        # minimizer_kwargs initialized to empty dictionary
        if len(minimizer_kwargs) == 0:
            self.minimizer_kwargs = {'method': 'Nelder-Mead',
                                     'options': {'disp': True,
                                                 'ftol': 1.0e-2,
                                                 'xtol': 1.0e-2}}
        else:

standard error!
    :param quantum_resource: quantum abstract machine object
    :param commutation_check: Optional flag toggling a safety check
                              ensuring all terms in `pauli_terms`
                              commute with each other
    :return: estimated expected value, expected value of each Pauli term in
             the sum, covariance matrix, variance of the estimator, and the
             number of shots taken.  The objected returned is a named tuple with
             field names as follows: expected_value, pauli_expectations,
             covariance, variance, n_shots.
             `expected_value' == coef_vec.dot(pauli_expectations)
    """
    if not isinstance(pauli_terms, (list, PauliSum)):
        raise TypeError("pauli_terms needs to be a list or a PauliSum")

    if isinstance(pauli_terms, PauliSum):
        pauli_terms = pauli_terms.terms

    # check if each term commutes with everything
    if commutation_check:
        if len(commuting_sets(sum(pauli_terms))) != 1:
            raise CommutationError("Not all terms commute in the expected way")

    program = copy.deepcopy(program)
    pauli_for_rotations = PauliTerm.from_list(
        [(value, key) for key, value in basis_transform_dict.items()])

    program += get_rotation_program(pauli_for_rotations)

    qubits = sorted(list(basis_transform_dict.keys()))
    for qubit in qubits:
        program.inst(MEASURE(qubit, qubit))

else:
            coeff = sum(t.coefficient for t in term_list)
            for t in term_list:
                if list(t._ops.items()) != list(first_term._ops.items()):
                    warnings.warn(
                        "The term {} will be combined with {}, but they have different "
                        "orders of operations. This doesn't matter for QVM or "
                        "wavefunction simulation but may be important when "
                        "running on an actual device.".format(
                            t.id(sort_ops=False), first_term.id(sort_ops=False)
                        )
                    )

            if not np.isclose(coeff, 0.0):
                terms.append(term_with_coeff(term_list[0], coeff))
    return PauliSum(terms)

:param program: pyquil Program preparing state
    :param quantum_resource: quantum abstract machine connection object
    :return: results dictionary where the key is the Pauli term ID and the value
             is the expected value
    """
    # group them into commuting sets and then measure
    grouped_terms = commuting_sets_by_zbasis(sum(pauli_terms))

    # measure the terms
    result_dictionary = {}
    for key, terms in grouped_terms.items():
        pauli_sum, identity_term = remove_identity(terms)
        if isinstance(identity_term, int):
            # no identity term
            pass
        elif isinstance(identity_term, PauliSum):
            result_dictionary[identity_term[0].id()] = 1.0
        else:
            print(identity_term, type(identity_term))
            raise TypeError("This type is not recognized for identity_term")

        results = estimate_pauli_sum(pauli_sum, dict(key), program,
                                     variance_bound / len(terms),
                                     quantum_resource)
        for idx, term in enumerate(pauli_sum.terms):
            result_dictionary[term.id()] = results.pauli_expectations[idx] / \
                                           term.coefficient
    return result_dictionary

:note: If the execution of ``quil_program`` is **non-deterministic**, i.e., if it includes
            measurements and/or noisy quantum gates, then the final wavefunction from which the
            expectation values are computed itself only represents a stochastically generated
            sample. The expectations returned from *different* ``expectation`` calls *will then
            generally be different*.

        :param Program prep_prog: Quil program for state preparation.
        :param Sequence[PauliTerm]|PauliSum pauli_terms: A list of PauliTerms or a PauliSum.
        :return: If ``pauli_terms`` is a PauliSum return its expectation value. Otherwise return
          a list of expectation values.
        :rtype: float|List[float]
        """

        is_pauli_sum = False
        if isinstance(pauli_terms, PauliSum):
            progs, coeffs = pauli_terms.get_programs()
            is_pauli_sum = True
        else:
            coeffs = [pt.coefficient for pt in pauli_terms]
            progs = [pt.program for pt in pauli_terms]

        bare_results = self.expectation(prep_prog, progs, needs_compilation=False, isa=False)
        results = [c * r for c, r in zip(coeffs, bare_results)]
        if is_pauli_sum:
            return sum(results)
        return results

def pauli_term_relabel(pauli_sum, label_map):
    """
    Relabel the elements of a pauli_sum via the `label_map`

    :param pauli_sum: pauli sum to relabel.  this can be a PauliTerm, PauliSum,
                      or even better a LIST!
    :param label_map: a dictionary mapping old label to new label
    :return: a list of pauli_terms relabeled
    """
    if isinstance(pauli_sum, PauliTerm):
        pauli_sum = PauliSum([pauli_sum])

    if isinstance(pauli_sum, PauliSum):
        pauli_sum = pauli_sum.terms

    relabeled_terms = []
    for term in pauli_sum:
        new_term_as_list = []
        for qlabel, pauli_element in term._ops.items():
            new_term_as_list.append((pauli_element, label_map[qlabel]))
        relabeled_terms.append(PauliTerm.from_list(
            new_term_as_list, coefficient=term.coefficient))
    return relabeled_terms

maxcut_graph.add_edge(*edge)
        graph = maxcut_graph.copy()

    cost_operators = []
    driver_operators = []
    if A is None or B is None:
        # bounds from arXiv: 1302:5843v3
        B = 1
        max_degree = max(graph.degree(), key=lambda x: graph.degree()[x])
        A_B_ratio = min(2*max_degree, len(graph.nodes()))/8.0
        A = A_B_ratio
        print A, B

    for node in graph.nodes():
        for node_j in graph.nodes():
            cost_operators.append(PauliSum([PauliTerm("Z", node, A)*PauliTerm("Z", node_j, 1.0)]))
        driver_operators.append(PauliSum([PauliTerm("X", node, -1)]))

    for i, j in graph.edges():
        cost_operators.append(PauliTerm("Z", i, -B*0.5)*PauliTerm("Z", j, 1.0) + PauliTerm("I", 0, B*0.5))

    # set number of qubits
    n_qubits = len(graph.nodes())

    if minimizer_kwargs is None:
        minimizer_kwargs = {'method': 'Nelder-Mead',
                            'options': {'ftol': 1.0e-2,
                                        'xtol': 1.0e-1,
                                        'disp': True,
                                        'maxfev': 100}}

    qaoa_inst = QAOA(CXN, n_qubits, steps=steps, cost_ham=cost_operators,

is then simplified according to the Pauli Algebra rules.

        :param other: a PauliSum, PauliTerm or Number object
        :return: A new PauliSum object given by the multiplication.
        :rtype: PauliSum
        """
        if not isinstance(other, (Number, PauliTerm, PauliSum)):
            raise ValueError(
                "Cannot multiply PauliSum by term that is not a Number, PauliTerm, or PauliSum"
            )
        elif isinstance(other, PauliSum):
            other_terms = other.terms
        else:
            other_terms = [other]
        new_terms = [lterm * rterm for lterm, rterm in product(self.terms, other_terms)]
        new_sum = PauliSum(new_terms)
        return new_sum.simplify()