Simulando circuitos classicamente com OBP

Neste guia, você vai aprender a simular classicamente instâncias de QuantumCircuit para estimar valores esperados inteiramente por meio da OBP.

Como a OBP pega um observável e o retropropaga por um circuito dado, a "simulação" de um circuito equivale a calcular o valor esperado do observável alvo em relação a esse circuito. Como você verá mais adiante, o pacote qiskit-addon-obp é capaz até de lidar com modelos de ruído simples, permitindo calcular também valores esperados com ruído!

Construindo um circuito de exemplo

Para os fins deste guia, usaremos o mesmo circuito de exemplo do guia de truncamento de termos de Pauli:

# Added by doQumentation — installs the packages this notebook needs if they are missing
import importlib.util

_needed = {"numpy": "numpy", "qiskit": "qiskit", "qiskit_addon_obp": "qiskit-addon-obp", "qiskit_addon_utils": "qiskit-addon-utils", "qiskit_ibm_runtime": "qiskit-ibm-runtime", "rustworkx": "rustworkx"}
_missing = [pip for module, pip in _needed.items()
            if importlib.util.find_spec(module) is None]
# One at a time, so a package that fails to install does not block the others
for _pip in _missing:
    %pip install -q {_pip}
if not _missing:
    print("\u2713 All required packages are installed")

import rustworkx.generators
from qiskit.synthesis import LieTrotter
from qiskit_addon_utils.problem_generators import (
    PauliOrderStrategy,
    generate_time_evolution_circuit,
    generate_xyz_hamiltonian,
)
from qiskit_addon_utils.slicing import combine_slices, slice_by_gate_types

# we generate a linear chain of 10 qubits
num_qubits = 10
linear_chain = rustworkx.generators.path_graph(num_qubits)

# we use an arbitrary XY model
hamiltonian = generate_xyz_hamiltonian(
    linear_chain,
    coupling_constants=(0.05, 0.02, 0.0),
    ext_magnetic_field=(0.02, 0.08, 0.0),
    pauli_order_strategy=PauliOrderStrategy.InteractionThenColor,
)
# we evolve for some time
circuit = generate_time_evolution_circuit(
    hamiltonian, synthesis=LieTrotter(reps=3), time=2.0
)
# slice the circuit by gate type
slices = slice_by_gate_types(circuit)

No entanto, o acima é puramente o circuito que descreve a evolução temporal sob um Hamiltoniano escolhido. Também precisamos de um estado inicial de onde partir, em relação ao qual calculamos os valores esperados do nosso observável.

Claro, poderíamos escolher o estado todo-zero (ou vácuo) como estado inicial, mas, para mostrar como inserir o próprio estado inicial, escolhemos um diferente abaixo.

Uma possibilidade seria preceder o circuito de evolução temporal acima com o estado inicial: circuit.compose(initial_state, front=True). Mas, como já dividimos o circuit em fatias, é mais fácil simplesmente inserir o estado inicial como a primeira fatia, o que fazemos abaixo.

Dessa forma, podemos simplesmente substituir a primeira fatia por outro estado inicial, se quisermos trocá-lo no futuro, sem precisar recalcular as fatias.

from qiskit.circuit import QuantumCircuit

initial_state = QuantumCircuit(num_qubits)
for i in range(0, num_qubits, 2):
    initial_state.x(i)

slices.insert(0, initial_state)

# for visualization purposes only, we recombine the slices with barriers between them and draw the resulting circuit
combine_slices(slices, include_barriers=True).draw("mpl", fold=50, scale=0.6)

Simulando um valor esperado sem ruído

Como observável alvo, escolhemos o observável ZZ nos qubits centrais:

from qiskit.quantum_info import SparsePauliOp

obs = SparsePauliOp("IIIIZZIIII")

Neste ponto, já estamos prontos para simular classicamente o valor esperado usando a OBP. Para isso, basta fornecer todas as fatias ao método backpropagate, assim:

from qiskit_addon_obp import backpropagate

vacuum_state_obs, _, metadata = backpropagate(obs, slices)

Agora retropropagamos nosso observável alvo obs por todo o circuito (incluindo o initial_state que colocamos em slices[0]), resultando em um novo SparsePauliOp cujo valor esperado obtemos projetando-o no estado de vácuo (|00...00>).

Isso pode ser feito de forma direta somando os coeficientes de todos os termos de Pauli definidos na base computacional:

vacuum_state_obs.coeffs[~vacuum_state_obs.paulis.x.any(axis=1)].sum()

np.complex128(-0.8285688012239535+4.9487770271457865e-20j)

Como verificação de sanidade (e para provar que isso funciona), podemos comparar nosso resultado com o Statevector do Qiskit:

from qiskit.quantum_info import Statevector

Statevector(combine_slices(slices)).expectation_value(obs)

np.complex128(-0.8285687255430366+0j)

Algumas observações sobre desempenho

A eficiência computacional da chamada de backpropagate acima dependerá fortemente de vários fatores, incluindo:

• a estrutura do circuit
• o método de divisão do circuito em fatias
• o observável alvo
• os parâmetros de truncamento

Como o método backpropagate simplifica o observável depois que cada fatia é aplicada, o número de portas em uma fatia pode influenciar dramaticamente a carga computacional. A estratégia mais agressiva em termos de simplificação do operador pode ser obtida dividindo seu circuito em fatias de portas individuais.

Além disso, você pode aproveitar todo o mecanismo de truncamento incorporado ao método backpropagate. Não fizemos isso acima, o que resulta efetivamente em um valor esperado exato, mas você pode aprender como fazê-lo no guia de truncamento de termos de Pauli.

Simulando um valor esperado com ruído

O pacote qiskit-addon-obp também oferece suporte ao tratamento de modelos de ruído na forma de PauliLindbladErrors. Isso é especialmente útil quando você caracterizou o modelo de ruído das camadas de 2 qubits do seu circuito, por exemplo usando o NoiseLearner.

Nesta seção, você verá como usar os objetos LayerError retornados pelo NoiseLearner para calcular valores esperados com ruído usando a OBP.

Obtendo um modelo de ruído

Normalmente, você executaria o NoiseLearner para obter um modelo de ruído do seu circuito específico. Para evitar complexidade (e aleatoriedade) neste tutorial, não faremos isso e, em vez disso, definiremos manualmente um modelo de ruído para o nosso circuito abaixo.

No entanto, garantimos que a estrutura dos nossos dados corresponda à do NoiseLearnerResult. Em sua forma atual (não transpilada), nosso circuito contém 4 camadas únicas de portas de 2 qubits:

• slices[1]: que tem portas Rxx atuando em todos os pares ímpares de qubits
• slices[2]: que tem portas Rxx atuando em todos os pares pares de qubits
• slices[3]: que tem portas Ryy atuando em todos os pares ímpares de qubits
• slices[4]: que tem portas Ryy atuando em todos os pares pares de qubits

Na célula abaixo, construímos manualmente 4 instâncias de LayerError, uma para cada uma dessas camadas, com algumas taxas de erro aleatórias.

import numpy as np
from qiskit.quantum_info import PauliList
from qiskit_ibm_runtime.utils.noise_learner_result import (
    LayerError,
    PauliLindbladError,
)

# fmt: off
pauli_errors_even = ['IIIIIIIIIX', 'IIIIIIIIIY', 'IIIIIIIIIZ', 'IIIIIIIIXI', 'IIIIIIIIXX', 'IIIIIIIIXY', 'IIIIIIIIXZ', 'IIIIIIIIYI', 'IIIIIIIIYX', 'IIIIIIIIYY', 'IIIIIIIIYZ', 'IIIIIIIIZI', 'IIIIIIIIZX', 'IIIIIIIIZY', 'IIIIIIIIZZ', 'IIIIIIIXII', 'IIIIIIIXXI', 'IIIIIIIXYI', 'IIIIIIIXZI', 'IIIIIIIYII', 'IIIIIIIYXI', 'IIIIIIIYYI', 'IIIIIIIYZI', 'IIIIIIIZII', 'IIIIIIIZXI', 'IIIIIIIZYI', 'IIIIIIIZZI', 'IIIIIIXIII', 'IIIIIIXXII', 'IIIIIIXYII', 'IIIIIIXZII', 'IIIIIIYIII', 'IIIIIIYXII', 'IIIIIIYYII', 'IIIIIIYZII', 'IIIIIIZIII', 'IIIIIIZXII', 'IIIIIIZYII', 'IIIIIIZZII', 'IIIIIXIIII', 'IIIIIXXIII', 'IIIIIXYIII', 'IIIIIXZIII', 'IIIIIYIIII', 'IIIIIYXIII', 'IIIIIYYIII', 'IIIIIYZIII', 'IIIIIZIIII', 'IIIIIZXIII', 'IIIIIZYIII', 'IIIIIZZIII', 'IIIIXIIIII', 'IIIIXXIIII', 'IIIIXYIIII', 'IIIIXZIIII', 'IIIIYIIIII', 'IIIIYXIIII', 'IIIIYYIIII', 'IIIIYZIIII', 'IIIIZIIIII', 'IIIIZXIIII', 'IIIIZYIIII', 'IIIIZZIIII', 'IIIXIIIIII', 'IIIXXIIIII', 'IIIXYIIIII', 'IIIXZIIIII', 'IIIYIIIIII', 'IIIYXIIIII', 'IIIYYIIIII', 'IIIYZIIIII', 'IIIZIIIIII', 'IIIZXIIIII', 'IIIZYIIIII', 'IIIZZIIIII', 'IIXIIIIIII', 'IIXXIIIIII', 'IIXYIIIIII', 'IIXZIIIIII', 'IIYIIIIIII', 'IIYXIIIIII', 'IIYYIIIIII', 'IIYZIIIIII', 'IIZIIIIIII', 'IIZXIIIIII', 'IIZYIIIIII', 'IIZZIIIIII', 'IXIIIIIIII', 'IXXIIIIIII', 'IXYIIIIIII', 'IXZIIIIIII', 'IYIIIIIIII', 'IYXIIIIIII', 'IYYIIIIIII', 'IYZIIIIIII', 'IZIIIIIIII', 'IZXIIIIIII', 'IZYIIIIIII', 'IZZIIIIIII', 'XIIIIIIIII', 'XXIIIIIIII', 'XYIIIIIIII', 'XZIIIIIIII', 'YIIIIIIIII', 'YXIIIIIIII', 'YYIIIIIIII', 'YZIIIIIIII', 'ZIIIIIIIII', 'ZXIIIIIIII', 'ZYIIIIIIII', 'ZZIIIIIIII']
pauli_errors_odd = ['IIIIIIIIXI', 'IIIIIIIIYI', 'IIIIIIIIZI', 'IIIIIIIXII', 'IIIIIIIXXI', 'IIIIIIIXYI', 'IIIIIIIXZI', 'IIIIIIIYII', 'IIIIIIIYXI', 'IIIIIIIYYI', 'IIIIIIIYZI', 'IIIIIIIZII', 'IIIIIIIZXI', 'IIIIIIIZYI', 'IIIIIIIZZI', 'IIIIIIXIII', 'IIIIIIXXII', 'IIIIIIXYII', 'IIIIIIXZII', 'IIIIIIYIII', 'IIIIIIYXII', 'IIIIIIYYII', 'IIIIIIYZII', 'IIIIIIZIII', 'IIIIIIZXII', 'IIIIIIZYII', 'IIIIIIZZII', 'IIIIIXIIII', 'IIIIIXXIII', 'IIIIIXYIII', 'IIIIIXZIII', 'IIIIIYIIII', 'IIIIIYXIII', 'IIIIIYYIII', 'IIIIIYZIII', 'IIIIIZIIII', 'IIIIIZXIII', 'IIIIIZYIII', 'IIIIIZZIII', 'IIIIXIIIII', 'IIIIXXIIII', 'IIIIXYIIII', 'IIIIXZIIII', 'IIIIYIIIII', 'IIIIYXIIII', 'IIIIYYIIII', 'IIIIYZIIII', 'IIIIZIIIII', 'IIIIZXIIII', 'IIIIZYIIII', 'IIIIZZIIII', 'IIIXIIIIII', 'IIIXXIIIII', 'IIIXYIIIII', 'IIIXZIIIII', 'IIIYIIIIII', 'IIIYXIIIII', 'IIIYYIIIII', 'IIIYZIIIII', 'IIIZIIIIII', 'IIIZXIIIII', 'IIIZYIIIII', 'IIIZZIIIII', 'IIXIIIIIII', 'IIXXIIIIII', 'IIXYIIIIII', 'IIXZIIIIII', 'IIYIIIIIII', 'IIYXIIIIII', 'IIYYIIIIII', 'IIYZIIIIII', 'IIZIIIIIII', 'IIZXIIIIII', 'IIZYIIIIII', 'IIZZIIIIII', 'IXIIIIIIII', 'IXXIIIIIII', 'IXYIIIIIII', 'IXZIIIIIII', 'IYIIIIIIII', 'IYXIIIIIII', 'IYYIIIIIII', 'IYZIIIIIII', 'IZIIIIIIII', 'IZXIIIIIII', 'IZYIIIIIII', 'IZZIIIIIII']
# fmt: on

np.random.seed(42)

layer_error_odd_xx = LayerError(
    circuit=slices[1],
    qubits=list(range(num_qubits)),
    error=PauliLindbladError(
        PauliList(pauli_errors_odd),
        0.0001 + 0.0004 * np.random.rand(len(pauli_errors_odd)),
    ),
)

layer_error_even_xx = LayerError(
    circuit=slices[2],
    qubits=list(range(num_qubits)),
    error=PauliLindbladError(
        PauliList(pauli_errors_even),
        0.0001 + 0.0004 * np.random.rand(len(pauli_errors_even)),
    ),
)

layer_error_odd_yy = LayerError(
    circuit=slices[3],
    qubits=list(range(num_qubits)),
    error=PauliLindbladError(
        PauliList(pauli_errors_odd),
        0.0001 + 0.0004 * np.random.rand(len(pauli_errors_odd)),
    ),
)

layer_error_even_yy = LayerError(
    circuit=slices[4],
    qubits=list(range(num_qubits)),
    error=PauliLindbladError(
        PauliList(pauli_errors_even),
        0.0001 + 0.0004 * np.random.rand(len(pauli_errors_even)),
    ),
)

Se você tivesse usado o NoiseLearner para identificar as camadas únicas de portas de 2 qubits do seu circuito e caracterizar o ruído delas, obteria um objeto NoiseLearnerResult. Esse resultado conteria uma lista de objetos LayerError, assim como os que construímos manualmente acima.

Para cada camada única de 2 qubits, o LayerError contém:

• o QuantumCircuit que representa essa camada de portas de 2 qubits
• os índices dos qubits sobre os quais esse circuito atua
• o PauliLindbladError que representa o modelo de ruído caracterizado dessa camada

O PauliLindbladError conterá dois objetos:

• a lista de erros de Pauli que foram caracterizados
• as taxas de erro correspondentes a cada um desses erros de Pauli

Normalmente, a lista de erros de Pauli será esparsa. Mais especificamente, ela conterá os erros de Pauli de um qubit em todos os qubits que têm portas atuando sobre eles, bem como os erros de Pauli de dois qubits em todos os qubits que estão conectados.

Inserindo as camadas ruidosas no nosso circuito

Na seção anterior, construímos especificamente um LayerError para cada uma das nossas camadas únicas de portas de 2 qubits conhecidas. Isso significa que sabemos exatamente qual LayerError corresponde a cada uma das nossas fatias.

Normalmente, ao usar o LayerError, você precisará descobrir qual é a camada única de portas de 2 qubits e onde ela ocorre dentro do seu circuito. Em seguida, será preciso ajustar seu circuit e/ou slices para inserir o LayerError de acordo. Como fazer isso no caso geral está fora do escopo deste guia prático. TODO: link to external documentation, once it exists! Aqui, nossa vida é mais simples porque sabemos a qual fatia um LayerError corresponde. Portanto, basta inserir novas fatias para representar o ruído.

Observe que precisamos envolver cada PauliLindbladError de LayerError.error em um PauliLindbladErrorInstruction para que seja uma instrução válida de QuantumCircuit.

from qiskit_addon_obp.utils.noise import PauliLindbladErrorInstruction

noisy_slices = []
for slice_ in slices:
    if slice_ == layer_error_even_xx.circuit:
        noisy_slices.append(
            QuantumCircuit.from_instructions(
                [
                    (
                        PauliLindbladErrorInstruction(
                            layer_error_even_xx.error
                        ),
                        slice_.qubits,
                    )
                ]
            )
        )
    elif slice_ == layer_error_odd_xx.circuit:
        noisy_slices.append(
            QuantumCircuit.from_instructions(
                [
                    (
                        PauliLindbladErrorInstruction(
                            layer_error_odd_xx.error
                        ),
                        slice_.qubits,
                    )
                ]
            )
        )
    elif slice_ == layer_error_even_yy.circuit:
        noisy_slices.append(
            QuantumCircuit.from_instructions(
                [
                    (
                        PauliLindbladErrorInstruction(
                            layer_error_even_yy.error
                        ),
                        slice_.qubits,
                    )
                ]
            )
        )
    elif slice_ == layer_error_odd_yy.circuit:
        noisy_slices.append(
            QuantumCircuit.from_instructions(
                [
                    (
                        PauliLindbladErrorInstruction(
                            layer_error_odd_yy.error
                        ),
                        slice_.qubits,
                    )
                ]
            )
        )
    noisy_slices.append(slice_)

Podemos conferir nosso trabalho e desenhar o circuito abaixo:

combine_slices(noisy_slices, include_barriers=True).draw(
    "mpl", fold=100, scale=0.8
)

Simulando um valor esperado com ruído

Neste ponto, simular classicamente o valor esperado funciona exatamente como antes, apenas

vacuum_state_noisy_obs, _, metadata = backpropagate(obs, noisy_slices)

vacuum_state_noisy_obs.coeffs[
    ~vacuum_state_noisy_obs.paulis.x.any(axis=1)
].sum()

np.complex128(-0.7230801696448901+7.082755280463563e-19j)

Destacamos novamente que várias questões de desempenho devem ser consideradas. Volte à seção correspondente acima.