Simulação do Hamiltoniano de Ising com circuitos dinâmicos

Estimativa de uso: 7.5 minutos em um processador Heron r3. (NOTA: Esta é apenas uma estimativa. O tempo de execução pode variar.)

Os circuitos dinâmicos são circuitos com feedforward clássico — em outras palavras, são medições no meio do circuito seguidas por operações lógicas clássicas que determinam operações quânticas condicionadas à saída clássica. Neste tutorial, simulamos o modelo Ising com efeito kick em uma rede hexagonal de spins e usamos circuitos dinâmicos para realizar interações além da conectividade física do hardware.

O modelo de Ising tem sido amplamente estudado em diversas áreas da física. Ele modela spins que sofrem interações de Ising entre pontos da rede, bem como impulsos do campo magnético local em cada ponto. A evolução temporal trotterizada dos spins considerados neste tutorial, retirada de [1], é dada pela seguinte unidade:

U(\theta)=\left(\prod_{\langle j, k\rangle} \exp \left(i \frac{\pi}{8} Z_j Z_k\right)\right)\left(\prod_j \exp \left(-i \frac{\theta}{2} X_j\right)\right)

Para investigar a dinâmica do spin, estudamos a magnetização média dos spins em cada local como uma função dos passos de Trotter. Portanto, construímos o seguinte observável:

\langle O\rangle = \frac{1}{N} \sum_i \langle Z_i \rangle

Para realizar a interação ZZ entre os pontos da rede, apresentamos uma solução utilizando o recurso de circuito dinâmico, levando a uma profundidade de dois qubits significativamente menor em comparação com o método de roteamento padrão com portas SWAP. Por outro lado, as operações clássicas de feedforward em circuitos dinâmicos normalmente têm tempos de execução mais longos do que os portões quânticos; portanto, os circuitos dinâmicos têm limitações e desvantagens. Também apresentamos uma maneira de adicionar uma sequência de desacoplamento dinâmico em qubits ociosos durante a operação clássica de feedforward usando a duração do alongamento.

Requisitos

Antes de iniciar este tutorial, certifique-se de que você tem o seguinte instalado:

Qiskit SDK v2.0 ou posterior com suporte para visualização
Qiskit Runtime v0.37 ou posterior com suporte à visualização (pip install 'qiskit-ibm-runtime[visualization]')
Biblioteca de gráficos Rustworkx (pip install rustworkx)
Qiskit Aer (pip install qiskit-aer)

Configurar

import numpy as np
from typing import List
import rustworkx as rx
import matplotlib.pyplot as plt
from rustworkx.visualization import mpl_draw
from qiskit.circuit import (
    Parameter,
    QuantumCircuit,
    QuantumRegister,
    ClassicalRegister,
)
from qiskit.transpiler import CouplingMap
from qiskit.quantum_info import SparsePauliOp
from qiskit.circuit.classical import expr
from qiskit.transpiler.preset_passmanagers import (
    generate_preset_pass_manager,
)
from qiskit.transpiler import PassManager
from qiskit.circuit.library import RZGate, XGate
from qiskit.transpiler.passes import (
    ALAPScheduleAnalysis,
    PadDynamicalDecoupling,
)

from qiskit.transpiler.basepasses import TransformationPass
from qiskit.circuit.measure import Measure
from qiskit.transpiler.passes.utils.remove_final_measurements import (
    calc_final_ops,
)
from qiskit.circuit import Instruction

from qiskit.visualization import plot_circuit_layout
from qiskit.circuit.tools import pi_check

from qiskit_aer import AerSimulator
from qiskit_aer.primitives import SamplerV2 as Aer_Sampler

from qiskit_ibm_runtime import (
    QiskitRuntimeService,
    Batch,
    SamplerV2 as Sampler,
)
from qiskit.providers.exceptions import QiskitBackendNotFoundError
from qiskit_ibm_runtime.visualization import (
    draw_circuit_schedule_timing,
)

Passo 1: Mapeie entradas clássicas para um circuito quântico

Começamos por definir a rede a simular. Optamos por trabalhar com a rede hexagonal (também chamada de favo de mel), que é um grafo plano com nós de grau 3. Aqui, especificamos o tamanho da rede, os parâmetros relevantes do circuito de interesse na dinâmica trotterizada. Simulamos a evolução temporal trotterizada sob o modelo de Ising sob três valores diferentes de $\theta$ do campo magnético local.

hex_rows = 3  # specify lattice size
hex_cols = 5
depths = range(9)  # specify Trotter steps
zz_angle = np.pi / 8  # parameter for ZZ interaction
max_angle = np.pi / 2  # max theta angle
points = 3  # number of theta parameters

θ = Parameter("θ")
params = np.linspace(0, max_angle, points)

def make_hex_lattice(hex_rows=1, hex_cols=1):
    """Define hexagon lattice."""
    hex_cmap = CouplingMap.from_hexagonal_lattice(
        hex_rows, hex_cols, bidirectional=False
    )
    data = list(hex_cmap.physical_qubits)
    graph = hex_cmap.graph.to_undirected(multigraph=False)
    edge_colors = rx.graph_misra_gries_edge_color(graph)
    layer_edges = {color: [] for color in edge_colors.values()}
    for edge_index, color in edge_colors.items():
        layer_edges[color].append(graph.edge_list()[edge_index])
    return data, layer_edges, hex_cmap, graph

Vamos começar com um pequeno exemplo de teste:

hex_rows_test = 1
hex_cols_test = 2

data_test, layer_edges_test, hex_cmap_test, graph_test = make_hex_lattice(
    hex_rows=hex_rows_test, hex_cols=hex_cols_test
)

# display a small example for illustration
node_colors_test = ["lightblue"] * len(graph_test.node_indices())
pos = rx.graph_spring_layout(
    graph_test,
    k=5 / np.sqrt(len(graph_test.nodes())),
    repulsive_exponent=1,
    num_iter=150,
)
mpl_draw(graph_test, node_color=node_colors_test, pos=pos)

Output:

Usaremos o pequeno exemplo para ilustração e simulação. Abaixo, também construímos um exemplo grande para mostrar que o fluxo de trabalho pode ser estendido para tamanhos grandes.

data, layer_edges, hex_cmap, graph = make_hex_lattice(
    hex_rows=hex_rows, hex_cols=hex_cols
)
num_qubits = len(data)
print(f"num_qubits = {num_qubits}")

# display the honeycomb lattice to simulate
node_colors = ["lightblue"] * len(graph.node_indices())
pos = rx.graph_spring_layout(
    graph,
    k=5 / np.sqrt(num_qubits),
    repulsive_exponent=1,
    num_iter=150,
)
mpl_draw(graph, node_color=node_colors, pos=pos)
plt.show()

Output:

num_qubits = 46

Construa circuitos unitários

Com o tamanho do problema e os parâmetros especificados, estamos agora prontos para construir o circuito parametrizado que simula a evolução temporal trotterizada de $U(\theta)$ com diferentes passos de Trotter, especificados pelo depth argumento. O circuito que construímos tem camadas alternadas de Rxportas ( $\theta$ ) e Rzz portas. Os Rzz portões realizam as interações ZZ entre spins acoplados, que serão colocados entre cada ponto da rede especificado pelo layer_edges argumento.

def gen_hex_unitary(
    num_qubits=6,
    zz_angle=np.pi / 8,
    layer_edges=[
        [(0, 1), (2, 3), (4, 5)],
        [(1, 2), (3, 4), (5, 0)],
    ],
    θ=Parameter("θ"),
    depth=1,
    measure=False,
    final_rot=True,
):
    """Build unitary circuit."""
    circuit = QuantumCircuit(num_qubits)
    # Build trotter layers
    for _ in range(depth):
        for i in range(num_qubits):
            circuit.rx(θ, i)
        circuit.barrier()
        for coloring in layer_edges.keys():
            for e in layer_edges[coloring]:
                circuit.rzz(zz_angle, e[0], e[1])
        circuit.barrier()
    # Optional final rotation, set True to be consistent with Ref. [1]
    if final_rot:
        for i in range(num_qubits):
            circuit.rx(θ, i)
    if measure:
        circuit.measure_all()

    return circuit

Visualize o pequeno circuito de teste:

circ_unitary_test = gen_hex_unitary(
    num_qubits=len(data_test),
    layer_edges=layer_edges_test,
    θ=Parameter("θ"),
    depth=1,
    measure=True,
)
circ_unitary_test.draw(output="mpl", fold=-1)

Output:

Da mesma forma, construa os circuitos unitários do exemplo grande em diferentes etapas de Trotter e o observável para estimar o valor esperado.

circuits_unitary = []
for depth in depths:
    circ = gen_hex_unitary(
        num_qubits=num_qubits,
        layer_edges=layer_edges,
        θ=Parameter("θ"),
        depth=depth,
        measure=True,
    )
    circuits_unitary.append(circ)

observables_unitary = SparsePauliOp.from_sparse_list(
    [("Z", [i], 1 / num_qubits) for i in range(num_qubits)],
    num_qubits=num_qubits,
)

Construa uma implementação de circuito dinâmico

Esta seção demonstra a implementação do circuito dinâmico principal para simular a mesma evolução temporal trotterizada. Observe que a rede alveolar que queremos simular não corresponde à rede pesada dos qubits de hardware. Uma maneira simples de mapear o circuito para o hardware é introduzir uma série de operações SWAP para aproximar os qubits interativos uns dos outros, a fim de realizar a interação ZZ. Aqui destacamos uma abordagem alternativa usando circuitos dinâmicos como solução, que ilustra que podemos usar a combinação de computação quântica e clássica em tempo real dentro de um circuito no Qiskit para realizar interações além do vizinho mais próximo.

Na implementação do circuito dinâmico, a interação ZZ é efetivamente implementada usando qubits auxiliares, medição no meio do circuito e feedforward. Para entender isso, observe que as rotações ZZ aplicam um fator de fase $e^{i\theta}$ ao estado com base em sua paridade. Para dois qubits, os estados da base computacional são $|00\rangle$ , $|01\rangle$ , $|10\rangle$ e $|11\rangle$ . A porta de rotação ZZ aplica um fator de fase aos estados $|01\rangle$ e $|10\rangle$ cuja paridade (o número de uns no estado) é ímpar e deixa os estados de paridade par inalterados. A seguir, descrevemos como podemos implementar efetivamente interações ZZ em dois qubits usando circuitos dinâmicos.

Calcule a paridade em um qubit auxiliar: em vez de aplicar diretamente ZZ a dois qubits, introduzimos um terceiro qubit, o qubit auxiliar, para armazenar as informações de paridade dos dois qubits de dados. Nós entrelaçamos a ancilla com cada qubit de dados usando portas CX do qubit de dados para o qubit da ancilla.
Aplique uma rotação Z de qubit único ao qubit auxiliar: isso ocorre porque o qubit auxiliar possui as informações de paridade dos dois qubits de dados, o que efetivamente implementa a rotação ZZ nos qubits de dados.
Meça o qubit auxiliar na base X: esta é a etapa fundamental que colapsa o estado do qubit auxiliar, e o resultado da medição nos diz o que aconteceu:
- Medida 0: quando se observa um resultado 0, significa que aplicamos corretamente uma rotação $ZZ(\theta)$ aos nossos qubits de dados.
- Medida 1: quando um resultado 1 é observado, aplicamos o modelo de regressão linear simples ( $ZZ(\theta + \pi)$ ) em seu lugar.
Aplicar porta de correção ao medir 1: Se medimos 1, aplicamos portas Z aos qubits de dados para “corrigir” a fase extra $\pi$ .

O circuito resultante é o seguinte:

Quando adotamos essa abordagem para simular uma rede em forma de colmeia, o circuito resultante se encaixa perfeitamente no hardware com uma rede hexagonal pesada: todos os qubits de dados residem nos locais e degree-3 es da rede, que forma uma rede hexagonal. Cada par de qubits de dados compartilha um qubit auxiliar localizado em um local de degree-2. Abaixo, construímos a rede de qubits para a implementação do circuito dinâmico, introduzindo qubits auxiliares (mostrados nos círculos roxos mais escuros).

def make_lattice(hex_rows=1, hex_cols=1):
    """Define heavy-hex lattice and corresponding lists of data and ancilla nodes."""
    hex_cmap = CouplingMap.from_hexagonal_lattice(
        hex_rows, hex_cols, bidirectional=False
    )
    data = list(hex_cmap.physical_qubits)

    heavyhex_cmap = CouplingMap()
    for d in data:
        heavyhex_cmap.add_physical_qubit(d)

    # make coupling map
    a = len(data)
    for edge in hex_cmap.get_edges():
        heavyhex_cmap.add_physical_qubit(a)
        heavyhex_cmap.add_edge(edge[0], a)
        heavyhex_cmap.add_edge(edge[1], a)
        a += 1
    ancilla = list(range(len(data), a))
    qubits = data + ancilla

    # color edges
    graph = heavyhex_cmap.graph.to_undirected(multigraph=False)
    edge_colors = rx.graph_misra_gries_edge_color(graph)
    layer_edges = {color: [] for color in edge_colors.values()}
    for edge_index, color in edge_colors.items():
        layer_edges[color].append(graph.edge_list()[edge_index])

    # construct observable
    obs_hex = SparsePauliOp.from_sparse_list(
        [("Z", [i], 1 / len(data)) for i in data],
        num_qubits=len(qubits),
    )

    return (data, qubits, ancilla, layer_edges, heavyhex_cmap, graph, obs_hex)

Visualize a rede hexagonal pesada para qubits de dados e qubits auxiliares em pequena escala:

(data, qubits, ancilla, layer_edges, heavyhex_cmap, graph, obs_hex) = (
    make_lattice(hex_rows=hex_rows, hex_cols=hex_cols)
)

print(f"number of data qubits = {len(data)}")
print(f"number of ancilla qubits = {len(ancilla)}")

node_colors = []
for node in graph.node_indices():
    if node in ancilla:
        node_colors.append("purple")
    else:
        node_colors.append("lightblue")

pos = rx.graph_spring_layout(
    graph,
    k=1 / np.sqrt(len(qubits)),
    repulsive_exponent=2,
    num_iter=200,
)

# Visualize the graph, blue circles are data qubits and purple circles are ancillas
mpl_draw(graph, node_color=node_colors, pos=pos)
plt.show()

Output:

number of data qubits = 46
number of ancilla qubits = 60

Abaixo, construímos o circuito dinâmico para a evolução temporal trotterizada. As RZZ portas são substituídas pela implementação do circuito dinâmico usando as etapas descritas acima.

def gen_hex_dynamic(
    depth=1,
    zz_angle=np.pi / 8,
    θ=Parameter("θ"),
    hex_rows=1,
    hex_cols=1,
    measure=False,
    add_dd=True,
):
    """Build dynamic circuits."""
    (data, qubits, ancilla, layer_edges, heavyhex_cmap, graph, obs_hex) = (
        make_lattice(hex_rows=hex_rows, hex_cols=hex_cols)
    )
    # Initialize circuit
    qr = QuantumRegister(len(qubits), "qr")
    cr = ClassicalRegister(len(ancilla), "cr")
    circuit = QuantumCircuit(qr, cr)

    for k in range(depth):
        # Single-qubit Rx layer
        for d in data:
            circuit.rx(θ, d)
        circuit.barrier()

        # CX gates from data qubits to ancilla qubits
        for same_color_edges in layer_edges.values():
            for e in same_color_edges:
                circuit.cx(e[0], e[1])
        circuit.barrier()

        # Apply Rz rotation on ancilla qubits and rotate into X basis
        for a in ancilla:
            circuit.rz(zz_angle, a)
            circuit.h(a)
        # Add barrier to align terminal measurement
        circuit.barrier()

        # Measure ancilla qubits
        for i, a in enumerate(ancilla):
            circuit.measure(a, i)
        d2ros = {}
        a2ro = {}
        # Retrieve ancilla measurement outcomes
        for a in ancilla:
            a2ro[a] = cr[ancilla.index(a)]

        # For each data qubit, retrieve measurement outcomes of neighboring
        # ancilla qubits
        for d in data:
            ros = [a2ro[a] for a in heavyhex_cmap.neighbors(d)]
            d2ros[d] = ros

        # Build classical feedforward operations (optionally add DD on idling
        # data qubits)
        for d in data:
            if add_dd:
                circuit = add_stretch_dd(circuit, d, f"data_{d}_depth_{k}")

            # # XOR the neighboring readouts of the data qubit;
            # if True, apply Z to it
            ros = d2ros[d]
            parity = ros[0]
            for ro in ros[1:]:
                parity = expr.bit_xor(parity, ro)
            with circuit.if_test(expr.equal(parity, True)):
                circuit.z(d)

        # Reset the ancilla if its readout is 1
        for a in ancilla:
            with circuit.if_test(expr.equal(a2ro[a], True)):
                circuit.x(a)
        circuit.barrier()

    # Final single-qubit Rx layer to match the unitary circuits
    for d in data:
        circuit.rx(θ, d)

    if measure:
        circuit.measure_all()
    return circuit, obs_hex


def add_stretch_dd(qc, q, name):
    """Add XpXm DD sequence."""
    s = qc.add_stretch(name)
    qc.delay(s, q)
    qc.x(q)
    qc.delay(s, q)
    qc.delay(s, q)
    qc.rz(np.pi, q)
    qc.x(q)
    qc.rz(-np.pi, q)
    qc.delay(s, q)
    return qc

Desacoplamento dinâmico (DD) e suporte para `stretch` duração

Uma ressalva ao uso da implementação de circuito dinâmico para realizar a interação ZZ é que a medição no meio do circuito e as operações clássicas de feedforward normalmente levam mais tempo para serem executadas do que os portões quânticos. Para suprimir a decoerência do qubit durante o tempo de inatividade para que as operações clássicas ocorram, adicionamos uma sequência de desacoplamento dinâmico (DD) após a operação de medição nos qubits auxiliares e antes da operação Z condicional no qubit de dados, antes da if_test instrução.

A sequência DD é adicionada pela função add_stretch_dd(), que usa o stretch duração para determinar os intervalos de tempo entre os portões DD. A stretch duração é uma forma de especificar um intervalo de tempo extensível para a delay operação, de modo que a duração do atraso possa aumentar para preencher o tempo de inatividade do qubit. As variáveis de duração especificadas por stretch são resolvidas em tempo de compilação nas durações desejadas que satisfazem uma determinada restrição. Isso é muito útil quando o tempo das sequências DD é essencial para obter um bom desempenho na supressão de erros. Para obter mais detalhes sobre o stretch tipo, consulte a OpenQASM documentação. Atualmente, o suporte para o stretch tipo em Qiskit Runtime é experimental. Para obter detalhes sobre as restrições de uso, consulte a seção de limitações da stretch documentação.

Usando as funções definidas acima, construímos os circuitos de evolução temporal trotterizados, com e sem DD, e os observáveis correspondentes.

Começamos visualizando o circuito dinâmico de um pequeno exemplo:

hex_rows_test = 1
hex_cols_test = 1

(
    data_test,
    qubits_test,
    ancilla_test,
    layer_edges_test,
    heavyhex_cmap_test,
    graph_test,
    obs_hex_test,
) = make_lattice(hex_rows=hex_rows_test, hex_cols=hex_cols_test)

node_colors = []
for node in graph_test.node_indices():
    if node in ancilla_test:
        node_colors.append("purple")
    else:
        node_colors.append("lightblue")
pos = rx.graph_spring_layout(
    graph_test,
    k=5 / np.sqrt(len(qubits_test)),
    repulsive_exponent=2,
    num_iter=150,
)

# display a small example for illustration
node_colors_test = ["lightblue"] * len(graph_test.node_indices())
mpl_draw(graph_test, node_color=node_colors, pos=pos)

Output:

circuit_dynamic_test, obs_dynamic_test = gen_hex_dynamic(
    depth=1,
    θ=Parameter("θ"),
    hex_rows=hex_rows_test,
    hex_cols=hex_cols_test,
    measure=False,
    add_dd=False,
)
circuit_dynamic_test.draw("mpl", fold=-1)

Output:

circuit_dynamic_dd_test, _ = gen_hex_dynamic(
    depth=1,
    θ=Parameter("θ"),
    hex_rows=hex_rows_test,
    hex_cols=hex_cols_test,
    measure=False,
    add_dd=True,
)
circuit_dynamic_dd_test.draw("mpl", fold=-1)

Output:

Da mesma forma, construa os circuitos dinâmicos para o exemplo grande:

circuits_dynamic = []
circuits_dynamic_dd = []
observables_dynamic = []
for depth in depths:
    circuit, obs = gen_hex_dynamic(
        depth=depth,
        θ=Parameter("θ"),
        hex_rows=hex_rows,
        hex_cols=hex_cols,
        measure=True,
        add_dd=False,
    )
    circuits_dynamic.append(circuit)

    circuit_dd, _ = gen_hex_dynamic(
        depth=depth,
        θ=Parameter("θ"),
        hex_rows=hex_rows,
        hex_cols=hex_cols,
        measure=True,
        add_dd=True,
    )
    circuits_dynamic_dd.append(circuit_dd)
    observables_dynamic.append(obs)

Etapa 2: Otimizar o problema para execução em hardware

Agora estamos prontos para transpilá-lo para o hardware. Vamos transpilá-las para o hardware, tanto a implementação padrão unitária quanto a implementação do circuito dinâmico.

Para transpilá-lo para o hardware, primeiro instanciamos o backend. Se disponível, escolheremos um backend que suporte a MidCircuitMeasure instrução (measure_2).

service = QiskitRuntimeService()
try:
    backend = service.least_busy(
        operational=True,
        simulator=False,
        use_fractional_gates=True,
        filters=lambda b: "measure_2" in b.supported_instructions,
    )
except QiskitBackendNotFoundError:
    backend = service.least_busy(
        operational=True,
        simulator=False,
        use_fractional_gates=True,
    )

Transpilação para circuitos dinâmicos

Primeiro, transpilamos os circuitos dinâmicos, com e sem adicionar a sequência DD. Para garantir que usamos o mesmo conjunto de qubits físicos em todos os circuitos para obter resultados mais consistentes, primeiro transpilamos o circuito uma vez e, em seguida, usamos seu layout para todos os circuitos subsequentes, especificados pelo initial_layout gerenciador de passagem. Em seguida, construímos os blocos unificados primitivos (PUBs) como entrada primitiva do Sampler.

pm_temp = generate_preset_pass_manager(
    optimization_level=3,
    backend=backend,
)
isa_temp = pm_temp.run(circuits_dynamic[-1])
dynamic_layout = isa_temp.layout.initial_index_layout(filter_ancillas=True)

pm = generate_preset_pass_manager(
    optimization_level=3, backend=backend, initial_layout=dynamic_layout
)

dynamic_isa_circuits = [pm.run(circ) for circ in circuits_dynamic]
dynamic_pubs = [(circ, params) for circ in dynamic_isa_circuits]

dynamic_isa_circuits_dd = [pm.run(circ) for circ in circuits_dynamic_dd]
dynamic_pubs_dd = [(circ, params) for circ in dynamic_isa_circuits_dd]

Podemos visualizar o layout do qubit do circuito transposto abaixo. Os círculos pretos mostram os qubits de dados e os qubits auxiliares usados na implementação do circuito dinâmico.

def _heron_coords_r2():
    cord_map = np.array(
        [
            [
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                3,
                7,
                11,
                15,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                1,
                5,
                9,
                13,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                3,
                7,
                11,
                15,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                1,
                5,
                9,
                13,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                3,
                7,
                11,
                15,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                1,
                5,
                9,
                13,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
                3,
                7,
                11,
                15,
                0,
                1,
                2,
                3,
                4,
                5,
                6,
                7,
                8,
                9,
                10,
                11,
                12,
                13,
                14,
                15,
            ],
            -1
            * np.array([j for i in range(15) for j in [i] * [16, 4][i % 2]]),
        ],
        dtype=int,
    )

    hcords = []
    ycords = cord_map[0]
    xcords = cord_map[1]
    for i in range(156):
        hcords.append([xcords[i] + 1, np.abs(ycords[i]) + 1])

    return hcords

plot_circuit_layout(
    dynamic_isa_circuits_dd[8],
    backend,
    qubit_coordinates=_heron_coords_r2(),
    view="virtual",
)

Output:

Note

Se você receber erros sobre neato não encontrado em plot_circuit_layout(), certifique-se de que o graphviz pacote está instalado e disponível em seu PATH. Se for instalado em um local não padrão (por exemplo, usando homebrew em MacOS ), talvez seja necessário atualizar sua variável PATH de ambiente. Isso pode ser feito dentro deste caderno usando o seguinte:

import os
os.environ['PATH'] = f"path/to/neato{os.pathsep}{os.environ['PATH']}"

dynamic_isa_circuits[1].draw(fold=-1, output="mpl", idle_wires=False)

Output:

dynamic_isa_circuits_dd[1].draw(fold=-1, output="mpl", idle_wires=False)

Output:

Transpilar usando `MidCircuitMeasure`

MidCircuitMeasure é um complemento às operações de medição disponíveis, calibrado especificamente para realizar medições no meio do circuito. A MidCircuitMeasure instrução corresponde à measure_2 instrução suportada pelos backends. Observe que não measure_2 é compatível com todos os backends. Você pode usar service.backends(filters=lambda b: "measure_2" in b.supported_instructions) para encontrar back-ends compatíveis com essa funcionalidade. Aqui, mostramos como transpilá-lo de modo que as medições no meio do circuito definidas no próprio circuito sejam executadas por meio da MidCircuitMeasure operação, caso o backend ofereça suporte a ela.

Abaixo, imprimimos a duração da measure_2 instrução e da instrução measure padrão.

print(
    f'Mid-circuit measurement `measure_2` duration: '
    f'{backend.instruction_durations.get('measure_2',0) * backend.dt * 1e9/1e3} μs'
)
print(
    f'Terminal measurement `measure` duration: '
    f'{backend.instruction_durations.get('measure',0) * backend.dt *1e9/1e3} μs'
)

Output:

Mid-circuit measurement `measure_2` duration: 1.3800000000000003 μs
Terminal measurement `measure` duration: 2.1800000000000006 μs

"""Pass that replaces terminal measures in the middle of the circuit with
MidCircuitMeasure instructions."""


class ConvertToMidCircuitMeasure(TransformationPass):
    """This pass replaces terminal measures in the middle of the circuit with
    MidCircuitMeasure instructions.
    """

    def __init__(self, target):
        super().__init__()
        self.target = target

    def run(self, dag):
        """Run the pass on a dag."""
        mid_circ_measure = None
        for inst in self.target.instructions:
            if isinstance(inst[0], Instruction) and inst[0].name.startswith(
                "measure_"
            ):
                mid_circ_measure = inst[0]
                break
        if not mid_circ_measure:
            return dag

        final_measure_nodes = calc_final_ops(dag, {"measure"})
        for node in dag.op_nodes(Measure):
            if node not in final_measure_nodes:
                dag.substitute_node(node, mid_circ_measure, inplace=True)

        return dag


pm = PassManager(ConvertToMidCircuitMeasure(backend.target))

dynamic_isa_circuits_meas2 = [pm.run(circ) for circ in dynamic_isa_circuits]
dynamic_pubs_meas2 = [(circ, params) for circ in dynamic_isa_circuits_meas2]

dynamic_isa_circuits_dd_meas2 = [
    pm.run(circ) for circ in dynamic_isa_circuits_dd
]
dynamic_pubs_dd_meas2 = [
    (circ, params) for circ in dynamic_isa_circuits_dd_meas2
]

Transpilação para circuitos unitários

Para estabelecer uma comparação justa entre os circuitos dinâmicos e sua contraparte unitária, usamos o mesmo conjunto de qubits físicos usados nos circuitos dinâmicos para os qubits de dados como layout para transpilagem dos circuitos unitários.

init_layout = [
    dynamic_layout[ind] for ind in range(circuits_unitary[0].num_qubits)
]


pm = generate_preset_pass_manager(
    target=backend.target,
    initial_layout=init_layout,
    optimization_level=3,
)


def transpile_minimize(circ: QuantumCircuit, pm: PassManager, iterations=10):
    """Transpile circuits for specified number of iterations and return the one
    with smallest two-qubit gate depth"""
    circs = [pm.run(circ) for i in range(iterations)]
    circs_sorted = sorted(
        circs,
        key=lambda x: x.depth(lambda x: x.operation.num_qubits == 2),
    )
    return circs_sorted[0]


unitary_isa_circuits = []
for circ in circuits_unitary:
    circ_t = transpile_minimize(circ, pm, iterations=100)
    unitary_isa_circuits.append(circ_t)

unitary_pubs = [(circ, params) for circ in unitary_isa_circuits]

Visualizamos o layout dos qubits dos circuitos unitários transpilados. Os círculos pretos indicam os qubits físicos usados para transpilá-los circuitos unitários e seus índices correspondem aos índices dos qubits virtuais. Ao comparar isso com o layout traçado para os circuitos dinâmicos, podemos confirmar que os circuitos unitários usam o mesmo conjunto de qubits físicos que os qubits de dados nos circuitos dinâmicos.

plot_circuit_layout(
    unitary_isa_circuits[-1],
    backend,
    qubit_coordinates=_heron_coords_r2(),
    view="virtual",
)

Output:

Agora adicionamos a sequência DD aos circuitos transpilados e construímos os PUBs correspondentes para o envio de tarefas.

pm_dd = PassManager(
    [
        ALAPScheduleAnalysis(target=backend.target),
        PadDynamicalDecoupling(
            dd_sequence=[
                XGate(),
                RZGate(np.pi),
                XGate(),
                RZGate(-np.pi),
            ],
            spacing=[1 / 4, 1 / 2, 0, 0, 1 / 4],
            target=backend.target,
        ),
    ]
)

unitary_isa_circuits_dd = pm_dd.run(unitary_isa_circuits)
unitary_pubs_dd = [(circ, params) for circ in unitary_isa_circuits_dd]

Compare a profundidade do portão de dois qubits dos circuitos unitários e dinâmicos

# compare circuit depth of unitary and dynamic circuit implementations
unitary_depth = [
    unitary_isa_circuits[i].depth(lambda x: x.operation.num_qubits == 2)
    for i in range(len(unitary_isa_circuits))
]

dynamic_depth = [
    dynamic_isa_circuits[i].depth(lambda x: x.operation.num_qubits == 2)
    for i in range(len(dynamic_isa_circuits))
]

plt.plot(
    list(range(len(unitary_depth))),
    unitary_depth,
    label="unitary circuits",
    color="#be95ff",
)
plt.plot(
    list(range(len(dynamic_depth))),
    dynamic_depth,
    label="dynamic circuits",
    color="#ff7eb6",
)
plt.xlabel("Trotter steps")
plt.ylabel("Two-qubit depth")
plt.legend()

Output:

<matplotlib.legend.Legend at 0x12628b0e0>

A principal vantagem do circuito baseado em medições é que, ao implementar múltiplas interações ZZ, as camadas CX podem ser paralelizadas e as medições podem ocorrer simultaneamente. Isso ocorre porque todas as interações ZZ comutam, portanto, o cálculo pode ser realizado com profundidade de medição 1. Após transpilar os circuitos, observamos que a abordagem do circuito dinâmico produz uma profundidade de dois qubits significativamente menor do que a abordagem unitária padrão, com a ressalva de que a medição adicional no meio do circuito e o feedforward clássico em si levam tempo e introduzem suas próprias fontes de erros.

Passo 3: Execute usando Qiskit primitives

Modo de teste local

Antes de enviar os trabalhos para o hardware, podemos executar uma pequena simulação de teste do circuito dinâmico usando o modo de teste local.

aer_sim = AerSimulator()
pm = generate_preset_pass_manager(backend=aer_sim, optimization_level=1)
circuit_dynamic_test.measure_all()
isa_qc = pm.run(circuit_dynamic_test)
with Batch(backend=aer_sim) as batch:
    sampler = Sampler(mode=batch)
    result = sampler.run([(isa_qc, params)]).result()

print(
    "Simulated average magnetization at trotter step = 1 at three theta values"
)
result[0].data["meas"].expectation_values(obs_dynamic_test[0])

Output:

Simulated average magnetization at trotter step = 1 at three theta values

array([ 0.16666667,  0.01529948, -0.14290365])

Simulação MPS

Para circuitos grandes, podemos usar o matrix_product_state simulador (MPS), que fornece um resultado aproximado do valor esperado de acordo com a dimensão da ligação escolhida. Posteriormente, utilizamos os resultados da simulação MPS como referência para comparar os resultados do hardware.

# The MPS simulation below took approximately 7 minutes to run on a
# laptop with Apple M1 chip

mps_backend = AerSimulator(
    method="matrix_product_state",
    matrix_product_state_truncation_threshold=1e-5,
    matrix_product_state_max_bond_dimension=100,
)
mps_sampler = Aer_Sampler.from_backend(mps_backend)

shots = 4096

data_sim = []
for j in range(points):
    circ_list = [
        circ.assign_parameters([params[j]]) for circ in circuits_unitary
    ]

    mps_job = mps_sampler.run(circ_list, shots=shots)
    result = mps_job.result()

    point_data = [
        result[d].data["meas"].expectation_values(observables_unitary)
        for d in depths
    ]

    data_sim.append(point_data)  # data at one theta value

data_sim = np.array(data_sim)

Com os circuitos e observáveis preparados, agora os executamos no hardware usando a primitiva Sampler.

Aqui, enviamos três trabalhos para unitary_pubs, dynamic_pubs, e dynamic_pubs_dd. Cada uma é uma lista de circuitos parametrizados correspondentes a nove etapas Trotter diferentes com três parâmetros diferentes de $\theta$ .

shots = 10000

with Batch(backend=backend) as batch:
    sampler = Sampler(mode=batch)

    sampler.options.experimental = {
        "execution": {
            "scheduler_timing": True
        },  # set to True to retrieve circuit timing info
    }

    job_unitary = sampler.run(unitary_pubs, shots=shots)
    print(f"unitary: {job_unitary.job_id()}")

    job_unitary_dd = sampler.run(unitary_pubs_dd, shots=shots)
    print(f"unitary_dd: {job_unitary_dd.job_id()}")

    job_dynamic = sampler.run(dynamic_pubs, shots=shots)
    print(f"dynamic: {job_dynamic.job_id()}")

    job_dynamic_dd = sampler.run(dynamic_pubs_dd, shots=shots)
    print(f"dynamic_dd: {job_dynamic_dd.job_id()}")

    job_dynamic_meas2 = sampler.run(dynamic_pubs_meas2, shots=shots)
    print(f"dynamic_meas2: {job_dynamic_meas2.job_id()}")

    job_dynamic_dd_meas2 = sampler.run(dynamic_pubs_dd_meas2, shots=shots)
    print(f"dynamic_dd_meas2: {job_dynamic_dd_meas2.job_id()}")

Output:

unitary: d96s4b52su3c739hakrg
unitary_dd: d96s4bt2su3c739haksg
dynamic: d96s4c0tcv6s73dk55mg
dynamic_dd: d96s4ckqp3as739qvid0
dynamic_meas2: d96s4csqp3as739qvie0
dynamic_dd_meas2: d96s4daf47jc73a5v8s0

Etapa 4: Pós-processamento e retorno dos resultados no formato clássico desejado

Após a conclusão das tarefas, podemos obter a duração do circuito a partir dos metadados dos resultados das tarefas e visualizar as informações da programação do circuito. Para saber mais sobre como visualizar as informações de programação de um circuito, consulte esta página.

# Circuit durations is reported in the unit of `dt`
# which can be retrieved from `Backend` object
unitary_durations = [
    job_unitary.result()[i].metadata["compilation"]["scheduler_timing"][
        "circuit_duration"
    ]
    for i in depths
]

dynamic_durations = [
    job_dynamic.result()[i].metadata["compilation"]["scheduler_timing"][
        "circuit_duration"
    ]
    for i in depths
]

dynamic_durations_meas2 = [
    job_dynamic_meas2.result()[i].metadata["compilation"]["scheduler_timing"][
        "circuit_duration"
    ]
    for i in depths
]

result_dd = job_dynamic_dd.result()[1]
circuit_schedule_dd = result_dd.metadata["compilation"]["scheduler_timing"][
    "timing"
]

# to visualize the circuit schedule, one can show the figure below
fig_dd = draw_circuit_schedule_timing(
    circuit_schedule=circuit_schedule_dd,
    included_channels=None,
    filter_readout_channels=False,
    filter_barriers=False,
    width=1000,
)

# Save to a file since the figure is large
fig_dd.write_html("scheduler_timing_dd.html")

Traçamos as durações dos circuitos unitários e dos circuitos dinâmicos. A partir do gráfico abaixo, podemos observar que, apesar do tempo necessário para medições no meio do circuito e operações clássicas, a implementação do circuito dinâmico com measure_2 resulta em durações de circuito comparáveis às da implementação unitária.

# visualize circuit durations


def convert_dt_to_microseconds(circ_duration: List, backend_dt: float):
    dt = backend_dt * 1e6  # dt in microseconds
    return list(map(lambda x: x * dt, circ_duration))


dt = backend.target.dt
plt.plot(
    depths,
    convert_dt_to_microseconds(unitary_durations, dt),
    color="#be95ff",
    linestyle=":",
    label="unitary",
)
plt.plot(
    depths,
    convert_dt_to_microseconds(dynamic_durations, dt),
    color="#ff7eb6",
    linestyle="-.",
    label="dynamic",
)
plt.plot(
    depths,
    convert_dt_to_microseconds(dynamic_durations_meas2, dt),
    color="#ff7eb6",
    linestyle="-.",
    marker="s",
    mfc="none",
    label="dynamic w/ meas2",
)

plt.xlabel("Trotter steps")
plt.ylabel(r"Circuit durations in $\mu$s")
plt.legend()

Output:

<matplotlib.legend.Legend at 0x12bfde270>

Após a conclusão das tarefas, recuperamos os dados abaixo e calculamos a magnetização média estimada pelas observáveis observables_unitary ou observables_dynamic que construímos anteriormente.

runs = {
    "unitary": (
        job_unitary,
        [observables_unitary] * len(circuits_unitary),
    ),
    "unitary_dd": (
        job_unitary_dd,
        [observables_unitary] * len(circuits_unitary),
    ),
    # Omitting Dyn w/o DD and Dynamic w/ DD plots for better readability
    # "dynamic": (job_dynamic, observables_dynamic),
    # "dynamic_dd": (job_dynamic_dd, observables_dynamic),
    "dynamic_meas2": (job_dynamic_meas2, observables_dynamic),
    "dynamic_dd_meas2": (
        job_dynamic_dd_meas2,
        observables_dynamic,
    ),
}

data_dict = {}
for key, (job, obs) in runs.items():
    data = []
    for i in range(points):
        data.append(
            [
                job.result()[ind].data["meas"].expectation_values(obs[ind])[i]
                for ind in depths
            ]
        )
    data_dict[key] = data

Abaixo, representamos graficamente a magnetização de spin em função dos passos de Trotter em diferentes valores de $\theta$ , correspondentes a diferentes intensidades do campo magnético local. Traçamos os resultados da simulação MPS pré-calculados para os circuitos ideais unitários, juntamente com os resultados experimentais a seguir:

operando os circuitos unitários com DD
operando os circuitos dinâmicos com DD e MidCircuitMeasure

plt.figure(figsize=(10, 6))

colors = ["#0f62fe", "#be95ff", "#ff7eb6"]
for i in range(points):
    plt.plot(
        depths,
        data_sim[i],
        color=colors[i],
        linestyle="solid",
        label=f"θ={pi_check(i*max_angle/(points-1))} (MPS)",
    )
    # plt.plot(
    #     depths,
    #     data_dict["unitary"][i],
    #     color=colors[i],
    #     linestyle=":",
    #     label=f"θ={pi_check(i*max_angle/(points-1))} (Unitary)",
    # )

    plt.plot(
        depths,
        data_dict["unitary_dd"][i],
        color=colors[i],
        marker="o",
        mfc="none",
        linestyle=":",
        label=f"θ={pi_check(i*max_angle/(points-1))} (Unitary w/DD)",
    )

    # Omitting Dyn w/o DD and Dynamic w/ DD plots for better readability
    # plt.plot(
    #     depths,
    #     data_dict["dynamic"][i],
    #     color=colors[i],
    #     linestyle="-.",
    #     label=f"θ={pi_check(i*max_angle/(points-1))} (Dyn w/o DD)",
    # )
    # plt.plot(
    #     depths,
    #     data_dict["dynamic_dd"][i],
    #     marker="D",
    #     mfc="none",
    #     color=colors[i],
    #     linestyle="-.",
    #     label=f"θ={pi_check(i*max_angle/(points-1))} (Dynamic w/ DD)",
    # )

    # plt.plot(
    #     depths,
    #     data_dict["dynamic_meas2"][i],
    #     color=colors[i],
    #     marker="s",
    #     mfc="none",
    #     linestyle=':',
    #     label=f"θ={pi_check(i*max_angle/(points-1))} (Dynamic w/ MidCircuitMeas)",
    # )

    plt.plot(
        depths,
        data_dict["dynamic_dd_meas2"][i],
        color=colors[i],
        marker="*",
        markersize=8,
        linestyle=":",
        label=f"θ={pi_check(i*max_angle/(points-1))} "
        f"(Dynamic w/ DD & MidCircuitMeas)",
    )


plt.xlabel("Trotter steps", fontsize=16)
plt.ylabel("Average magnetization", fontsize=16)
plt.xticks(rotation=45)
handles, labels = plt.gca().get_legend_handles_labels()
plt.legend(
    handles,
    labels,
    loc="upper right",
    bbox_to_anchor=(1.46, 1.0),
    shadow=True,
    ncol=1,
)
plt.title(
    f"{hex_rows}x{hex_cols} hex ring, {num_qubits} data qubits, "
    f"{len(ancilla)} ancilla qubits \n{backend.name}: Sampler"
)
plt.show()

Output:

Quando comparamos os resultados experimentais com a simulação, vemos que a implementação do circuito dinâmico (linha pontilhada com estrelas) tem, em geral, um desempenho melhor do que a implementação unitária padrão (linha pontilhada com círculos). Em resumo, apresentamos circuitos dinâmicos como uma solução para simular modelos de spin de Ising em uma rede hexagonal, uma topologia que não é nativa do hardware. A solução de circuito dinâmico permite interações ZZ entre qubits que não são vizinhos mais próximos, com uma profundidade de porta de dois qubits mais curta do que usando portas SWAP, ao custo de introduzir qubits auxiliares extras e operações clássicas de feedforward.

Referências

[1] Computação quântica com Qiskit, por Javadi-Abhari, A., Treinish, M., Krsulich, K., Wood, C.J Lishman, J., Gacon, J., Martiel, S., Nation, P.D Bispo, L.S Cross, A.W. e Johnson, B.R 2024. Pré-impressão de “ arXiv ” (A evolução da ciência da saúde em um mundo em mudança) arXiv:2405.08810 (2024)

Esta página foi útil?

Relate um bug, erro de digitação ou solicite conteúdo no GitHub.

Simulação do Hamiltoniano de Ising com circuitos dinâmicos

Requisitos

Configurar

Passo 1: Mapeie entradas clássicas para um circuito quântico

Construa circuitos unitários

Construa uma implementação de circuito dinâmico

Desacoplamento dinâmico (DD) e suporte para stretch duração

Etapa 2: Otimizar o problema para execução em hardware

Transpilação para circuitos dinâmicos

Transpilar usando MidCircuitMeasure

Transpilação para circuitos unitários

Compare a profundidade do portão de dois qubits dos circuitos unitários e dinâmicos

Passo 3: Execute usando Qiskit primitives

Modo de teste local

Simulação MPS

Etapa 4: Pós-processamento e retorno dos resultados no formato clássico desejado

Referências

Desacoplamento dinâmico (DD) e suporte para `stretch` duração

Transpilar usando `MidCircuitMeasure`