def METHOD_NAME(self, evolution):
    # get operators and time to evolve
    operators = evolution.operator
    time = evolution.time
    if not isinstance(operators, list):
        pauli_list = [(Pauli(op), coeff) for op, coeff in operators.to_list()]
        coeffs = [np.real(coeff) for op, coeff in operators.to_list()]
    else:
        pauli_list = [(op, 1) for op in operators]
        coeffs = [1 for op in operators]
    # We artificially make the weights positive
    weights = np.abs(coeffs)
    lambd = np.sum(weights)
    num_gates = int(np.ceil(2 * (lambd**2) * (time**2) * self.reps))
    # The protocol calls for the removal of the individual coefficients,
    # and multiplication by a constant evolution time.
    evolution_time = lambd * time / num_gates
    self.sampled_ops = algorithm_globals.random.choice(
        np.array(pauli_list, dtype=object),
        size=(num_gates,),
        p=weights / lambd,
    )
    # Update the coefficients of sampled_ops
    self.sampled_ops = [(op, evolution_time) for op, coeff in self.sampled_ops]
    # pylint: disable=cyclic-import
    from qiskit.circuit.library.pauli_evolution import PauliEvolutionGate
    # Build the evolution circuit using the LieTrotter synthesis with the sampled operators
    lie_trotter = LieTrotter(
        insert_barriers=self.insert_barriers, atomic_evolution=self.atomic_evolution
    )
    evolution_circuit = PauliEvolutionGate(
        sum(SparsePauliOp(op) for op, coeff in self.sampled_ops),
        time=evolution_time,
        synthesis=lie_trotter,
    ).definition
    return evolution_circuit