Fórmulas multiprodutos para reduzir o erro de Trotter

Estimativa de tempo de execução: Quatro minutos em um processador Heron r2 (NOTA: Trata-se apenas de uma estimativa. (O tempo de execução pode variar.)

Resultados do aprendizado

• Como as fórmulas multiproduto (MPFs) reduzem o erro de Trotter na simulação hamiltoniana ao combinar valores esperados de vários circuitos rasos
• Quando as MPFs são mais vantajosas do que as fórmulas padrão dos produtos e quando não são a ferramenta adequada
• Como calcular os coeficientes MPF estáticos e dinâmicos usando o qiskit_addon_mpf pacote
• Como executar um fluxo de trabalho do MPF de ponta a ponta em um hardwar IBM Quantum®, incluindo transpilação, mitigação de erros e pós-processamento

Pré-requisitos

• Métodos de compilação para circuitos de simulação hamiltoniana — apresenta os circuitos de Trotter (fórmula do produto) no Qiskit.
• Fórmulas de produtos no Qiskit, em especial as SuzukiTrotter classes de síntese e LieTrotter .
• Qiskit primitives e a interface do Estimador.

Segundo plano

O que são fórmulas multiprodutos?

Ao simular sistemas quânticos em um computador quântico, uma tarefa central consiste em aproximar o operador de evolução temporal $e^{-iHt}$ para um hamiltoniano $H$. A abordagem padrão utiliza fórmulas de produto (PFs), também conhecidas como decomposições de Trotter-Suzuki. Essas funções decompõem $H = \sum_{a=1}^d F_a$ em termos cujas operações unitárias individuais $e^{-iF_a t}$ são eficientes de implementar e, em seguida, aproximam a evolução completa como um produto ordenado dessas operações unitárias mais simples.

A fórmula do produto de primeira ordem (Lie-Trotter) é:

o que resulta em um erro quadrático: $S_1(t) = e^{-iHt} + \mathcal{O}(t^2)$. Fórmulas simétricas de ordem superior $S_{2\chi}(t)$, em que $\chi$ indica a ordem da fórmula do produto simétrico (ver Ref. [1] ), convergem mais rapidamente, conforme $e^{-iHt} + \mathcal{O}(t^{2\chi+1})$, mas à custa de circuitos mais profundos a cada passo.

Para reduzir o erro em uma ordem fixa $\chi$, costuma-se dividir o tempo total de evolução $t$ em $k$ passos de Trotter menores. Cada etapa aproxima $e^{-iHt/k}$ por meio de uma fórmula de produto, e as etapas são concatenadas:

Para uma fórmula simétrica de ordem $2\chi$, o erro residual de Trotter varia proporcionalmente a $\mathcal{O}\!\left(t^{2\chi+1} / k^{2\chi}\right)$. Assim, aumentar $k$ suprime rapidamente o erro de Trotter — mas também aumenta linearmente a profundidade do circuito, e em hardware sujeito a ruído isso significa mais ruído acumulado nas portas. Essa tensão entre o erro de Trotter (que favorece valores maiores de $k$ ) e o ruído do hardware (que favorece valores menores de $k$ ) é exatamente o que as fórmulas multiproduto se propõem a resolver. Observe que as MPFs consistem em combinar resultados de diferentes escolhas de $k$ em uma ordem fixa $\chi$ — elas não alteram a ordem da fórmula do produto subjacente.

As fórmulas multiproduto (MPFs) [1] constroem uma combinação linear ponderada dos valores esperados obtidos a partir de vários circuitos de Trotter menos profundos, cada um utilizando um número diferente de etapas de Trotter $k_1, k_2, \ldots, k_r$ (um conjunto de contagens de etapas $r$ ):

onde $\langle A \rangle_{k_j}(t)$ é o valor esperado de um observável $A$ no instante $t$, estimado a partir de um circuito de Trotter com $k_j$ passos, e os coeficientes $\{x_j\}_{j=1}^r$ são escolhidos de modo que os principais termos de erro de Trotter na combinação se cancelem. Voltaremos a essa expressão na Etapa 4, onde a avaliaremos explicitamente para combinar nossos resultados de Trotter. O ponto prático fundamental é que o circuito mais profundo do MPF precisa apenas de $k_{\max}$ etapas, o que é muito menor do que o único $k$ que seria necessário para atingir diretamente o mesmo erro efetivo de Trotter. Os circuitos menos complexos tornam a abordagem MPF mais adequada para hardware com ruído.

Como os coeficientes são determinados?

Existem duas famílias de coeficientes MPF:

Os coeficientes estáticos são independentes do hamiltoniano, do estado inicial e do tempo de evolução. Elas são determinadas pela resolução de um sistema linear $Ax = b$ que garante o cancelamento dos principais termos de erro de Trotter. Para um conjunto de passos de Trotter $\{k_j\}_{j=1}^r$ utilizado com uma fórmula de produto simétrico de ordem $2\chi$, a expansão do erro de Trotter em potências inversas de $k_j$ leva a equações de restrição da forma:

onde os expoentes inteiros $\{\eta_n\}$ são as ordens dos termos sucessivos do erro de Trotter para a fórmula de produto escolhida. Para um PF simétrico de ordem $2\chi$, o erro principal em $\left[S_{2\chi}(t/k)\right]^k$ varia proporcionalmente a $1/k^{2\chi}$, com correções subsequentes em $1/k^{2\chi+2}, 1/k^{2\chi+4}, \ldots$ — portanto, os expoentes são $\eta_n = 2\chi + 2n$. Para PFs não simétricos, tanto as potências ímpares quanto as pares contribuem e $\eta_n = 2\chi + n$. Consulte a Ref. [1] para a derivação completa. A primeira equação do sistema acima garante a imparcialidade (o MPF reproduz o valor exato da expectativa no limite $k_j \to \infty$ ), e as demais equações $r-1$ cancelam sucessivamente os primeiros termos de erro de Trotter $r-1$. Quando a norma $L_1$ resultante $\|x\|_1$ é muito grande (o que amplifica o ruído de amostragem), é possível, em vez disso, resolver um problema de otimização aproximada que limite $\|x\|_1$ e, ao mesmo tempo, minimize $\|Ax - b\|$.

Os coeficientes dinâmicos [2], [3] dependem, além disso, do hamiltoniano, do estado inicial e do tempo de evolução $t$. Eles minimizam a distância na norma de Frobenius entre o estado real evoluído no tempo e a aproximação MPF:

onde $M_{ij}(t) = \mathrm{Tr}[\rho_{k_i}(t)\,\rho_{k_j}(t)]$ é a matriz de Gram das sobreposições entre os estados evoluídos por Trotter para diferentes contagens de passos $k_i, k_j$, e $L_i(t) = \mathrm{Tr}[\rho(t)\,\rho_{k_i}(t)]$ mede a sobreposição com o estado exato (aproximado). Neste tutorial, essas grandezas são calculadas de forma eficiente utilizando métodos de redes tensoriais, especificamente os backends TeNPy-based em qiskit_addon_mpf.

Quando usar os MPFs

Os MPFs são mais benéficos quando:

• A profundidade do circuito é o gargalo. Se o ruído do hardware limitar a profundidade em que é possível operar, utilize MPFs para obter maior precisão efetiva do método de Trotter em circuitos menos profundos.
• Você precisa de valores de expectativa precisos, e não da preparação completa do estado. As MPFs operam no nível dos valores esperados — elas combinam números clássicos, e não estados quânticos. Portanto, são ideais para estimativa observável ao se utilizar a primitiva Estimador.
• Você combina um número modesto de passos de trot. Normalmente, a combinação de $r = 3$ – $5$, com diferentes contagens de passos $k_j$, é suficiente para eliminar vários termos de erro de Trotter de ordem superior, mantendo $\|x\|_1$ em um nível gerenciável.

Quando os MPFs podem não ser úteis

• Tempos de evolução muito curtos. Quando $t$ é pequeno o suficiente para que uma única fórmula de Trotter de baixa ordem já seja precisa, o esforço adicional de executar vários circuitos torna-se desnecessário.
• Tarefas de preparação para o exame estadual. Os MPFs produzem um valor esperado corrigido, e não um estado quântico corrigido. Se você precisar do estado real ao longo do tempo (por exemplo, como entrada para outra sub-rotina quântica), os MPFs não se aplicam.
• Contagens de passos do método Trotter que violam o regime de convergência. A derivação do coeficiente estático expande cada “ $\left[S_{2\chi}(t/k_j)\right]^{k_j}$ ” individual como uma série em “ $t/k_j$ ”; essa expansão só converge bem quando “ $t/k_{\min} \lesssim 1$ ”. Se “ $k_{\min}$ ” for escolhido muito pequeno para o “ $t$ ” dado, o circuito mais superficial fica bem fora do regime perturbativo, os termos de erro de ordem superior que o MPF deixa sem cancelar tornam-se grandes, e o cancelamento pode exigir coeficientes grandes. A norma de $L_1$ $\|x\|_1$ serve como um diagnóstico prático: quando $\|x\|_1 \gg 1$, a sobrecarga de amostragem $\propto \|x\|_1^2$ pode superar a redução do erro de Trotter. Consulte o guia sobre como escolher degraus Trotter para obter mais detalhes.

O que este tutorial aborda

Este tutorial apresenta um fluxo de trabalho completo do MPF em duas etapas. Primeiro, um exemplo de simulador em pequena escala (cadeia de Heisenberg de 10 qubits) demonstra como definir o problema, calcular os coeficientes MPF estáticos e dinâmicos e comparar os valores esperados resultantes com a diagonalização exata. Em seguida, um exemplo de hardware em grande escala (cadeia XXZ de 50 qubits) mostra como fazer a transpilagem, executar em um hardwar IBM Quantum e com mitigação de erros e processar os resultados posteriormente usando os coeficientes MPF. Ao longo do texto, utilizamos o qiskit_addon_mpf pacote em conjunto com as ferramentas padrão do Qiskit.

Requisitos

Antes de iniciar este tutorial, verifique se você tem os seguintes itens instalados:

• Qiskit SDK v2.0 ou versão posterior, com suporte à visualização
• Qiskit Runtime v0.22 ou posterior (pip install qiskit-ibm-runtime)
• Simulador Qiskit Aer (pip install qiskit-aer)
• Complemento MPF Qiskit com o backend TeNPy (pip install "qiskit-addon-mpf[tenpy]")
• Utilitários do complemento Qiskit (pip install qiskit-addon-utils)
• SciPy (pip install scipy)

Instalação

A seguir, reunimos todas as importações de pacotes utilizadas ao longo deste tutorial em uma única célula. Também definimos uma CollectAndCollapse etapa do transpiler que funde rotações adjacentes rxx e ryy em uma única XXPlusYYGaterotação. Essa etapa é aplicada tanto durante a construção do circuito na Etapa 1 (para manter baixo o número de portas) quanto indiretamente quando extraímos a estrutura em camadas para o MPF dinâmico na Etapa 4 (o TeNPy espera portas de dois qubits, e não pares de rotações não fundidas).

import warnings

import numpy as np
import matplotlib.pyplot as plt
from functools import partial
from copy import deepcopy

from qiskit import QuantumCircuit
from qiskit.quantum_info import Pauli, SparsePauliOp, Statevector
from qiskit.synthesis import SuzukiTrotter
from qiskit.transpiler import CouplingMap, PassManager
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
from qiskit.circuit.library import XXPlusYYGate
from qiskit.transpiler.passes.optimization.collect_and_collapse import (
    CollectAndCollapse,
    collect_using_filter_function,
    collapse_to_operation,
)

from qiskit_aer import AerSimulator
from qiskit_ibm_runtime import EstimatorV2 as Estimator, QiskitRuntimeService

from qiskit_addon_utils.problem_generators import (
    generate_xyz_hamiltonian,
    generate_time_evolution_circuit,
)
from qiskit_addon_utils.slicing import slice_by_depth
from qiskit_addon_mpf.static import setup_static_lse
from qiskit_addon_mpf.dynamic import setup_dynamic_lse
from qiskit_addon_mpf.costs import (
    setup_exact_problem,
    setup_sum_of_squares_problem,
    setup_frobenius_problem,
)
from qiskit_addon_mpf.backends.tenpy_layers import (
    LayerModel,
    LayerwiseEvolver,
)
from qiskit_addon_mpf.backends.tenpy_tebd import MPOState, MPS_neel_state

from scipy.linalg import expm

# Suppress TeNPy's `unit_cell_width` future-API warning. The default
# (`unit_cell_width=len(sites)`) is correct for Chain lattices, which is what
# `CouplingMap.from_line(...)` produces here, so the warning is informational.
warnings.filterwarnings(
    "ignore",
    message=r".*unit_cell_width.*",
    category=UserWarning,
)


# --- Helper: collect XX + YY rotations into a single gate ---
def filter_function(node):
    return node.op.name in {"rxx", "ryy"}


collect_function = partial(
    collect_using_filter_function,
    filter_function=filter_function,
    split_blocks=True,
    min_block_size=1,
)


def collapse_to_xx_plus_yy(block):
    param = 0.0
    for node in block.data:
        param += node.operation.params[0]
    return XXPlusYYGate(param)


collapse_function = partial(
    collapse_to_operation,
    collapse_function=collapse_to_xx_plus_yy,
)

pm = PassManager()
pm.append(CollectAndCollapse(collect_function, collapse_function))

Exemplo de simulador em pequena escala

Passo 1: Mapear entradas clássicas para um problema quântico

Começamos com um modelo de Heisenberg de 10 qubits em uma reta, utilizando o estado de Néel $\vert 0101\ldots01 \rangle$ como estado inicial. O hamiltoniano é:

onde $J$ é a intensidade do acoplamento entre vizinhos mais próximos. Medimos o correlacionador ZZ $Z_{L/2-1} Z_{L/2}$ em um par de qubits no meio da cadeia e utilizamos os passos de Trotter $k_j = [1, 2, 4]$ com uma fórmula de produto de segunda ordem.

L = 10

# Generate coupling map and Hamiltonian
coupling_map = CouplingMap.from_line(L, bidirectional=False)

hamiltonian = generate_xyz_hamiltonian(
    coupling_map,
    coupling_constants=(1.0, 1.0, 1.0),
    ext_magnetic_field=(0.0, 0.0, 0.0),
)
print(hamiltonian)

Output:

SparsePauliOp(['IIIIIIIXXI', 'IIIIIIIYYI', 'IIIIIIIZZI', 'IIIIIXXIII', 'IIIIIYYIII', 'IIIIIZZIII', 'IIIXXIIIII', 'IIIYYIIIII', 'IIIZZIIIII', 'IXXIIIIIII', 'IYYIIIIIII', 'IZZIIIIIII', 'IIIIIIIIXX', 'IIIIIIIIYY', 'IIIIIIIIZZ', 'IIIIIIXXII', 'IIIIIIYYII', 'IIIIIIZZII', 'IIIIXXIIII', 'IIIIYYIIII', 'IIIIZZIIII', 'IIXXIIIIII', 'IIYYIIIIII', 'IIZZIIIIII', 'XXIIIIIIII', 'YYIIIIIIII', 'ZZIIIIIIII'],
              coeffs=[1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j,
 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j,
 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j])

# Observable: ZZ on the middle pair of qubits
observable = SparsePauliOp.from_sparse_list(
    [("ZZ", (L // 2 - 1, L // 2), 1.0)], num_qubits=L
)
print(observable)

Output:

SparsePauliOp(['IIIIZZIIII'],
              coeffs=[1.+0.j])

# MPF parameters
mpf_trotter_steps = [1, 2, 4]
order = 2
symmetric = False

trotter_times = np.arange(0.5, 1.55, 0.1)
exact_evolution_times = np.arange(trotter_times[0], 1.55, 0.05)

Construir circuitos de Trotter

Criamos os circuitos implementando as evoluções temporais aproximadas de Trotter para cada ponto no tempo e cada número de passos de Trotter. A CollectAndCollapse etapa definida na seção “Configuração” agrupa as rotações XX e YY em portas únicas do tipo XX+YY, a fim de preparar uma simulação mais eficiente da rede tensorial posteriormente.

# Initial Neel state preparation
initial_state_circ = QuantumCircuit(L)
initial_state_circ.x([i for i in range(L) if i % 2 != 0])


all_circs = []
for total_time in trotter_times:
    mpf_trotter_circs = [
        generate_time_evolution_circuit(
            hamiltonian,
            time=total_time,
            synthesis=SuzukiTrotter(reps=num_steps, order=order),
        )
        for num_steps in mpf_trotter_steps
    ]

    mpf_trotter_circs = pm.run(
        mpf_trotter_circs
    )  # Collect XX and YY into XX + YY

    mpf_circuits = [
        initial_state_circ.compose(circuit) for circuit in mpf_trotter_circs
    ]
    all_circs.append(mpf_circuits)

mpf_circuits[-1].draw("mpl", fold=-1)

Output:

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

Para o exemplo em pequena escala, vamos usar o simulador Aer. Duas transformações ocorrem antes que os circuitos estejam prontos para serem executados:

1. Coleta de portas no nível da simulação hamiltoniana. XXPlusYYGateNa célula “Setup”, criamos uma CollectAndCollapse passagem que combina as rotações adjacentes rxx e ryy em uma única. Já aplicamos essa etapa quando construímos os circuitos de Trotter na Etapa 1 (a pm.run(...) chamada). Isso reduz o número de portas de dois qubits e, ao mesmo tempo, gera uma estrutura mais adequada para a simulação por rede tensorial no cálculo posterior dos coeficientes dinâmicos.

2. Conversão para o ISA do simulador. A seguir, executamos o gerenciador de passagens com a predefinição do Qiskit para optimization_level=3 adaptar cada circuito de Trotter à arquitetura do conjunto de instruções (ISA) do simulador.

aer_sim = AerSimulator()
pm_sim = generate_preset_pass_manager(backend=aer_sim, optimization_level=3)

isa_circs_all_times = [
    pm_sim.run([deepcopy(c) for c in mpf_circuits])
    for mpf_circuits in all_circs
]

Passo 3: Execute usando Qiskit primitives

Para o exemplo em pequena escala, executamos os circuitos de Trotter reduzidos por ISA por meio da EstimatorV2 primitiva implementada pelo Aer. Isso nos proporciona um valor de referência sem ruído para cada par $(k_j, t)$ — esses são os valores $\langle A \rangle_{k_j}(t)$ que o MPF combinará na Etapa 4. Analisamos os períodos de evolução para que, posteriormente, possamos traçar a curva completa da série temporal de cada fórmula de produto individual e da MPF.

estimator = Estimator(mode=aer_sim)

mpf_expvals_all_times, mpf_stds_all_times = [], []
for isa_circuits in isa_circs_all_times:
    result = estimator.run(
        [(circuit, observable) for circuit in isa_circuits], precision=0.005
    ).result()
    mpf_expvals_all_times.append([res.data.evs for res in result])
    mpf_stds_all_times.append([res.data.stds for res in result])

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

A Etapa 4 é onde o MPF é efetivamente construído. Embora os coeficientes $x_j$ sejam calculados aqui (e, para a variante dinâmica, esse cálculo possa ser intensivo), conceitualmente eles constituem uma fórmula clássica para combinar as medições quânticas da Etapa 3 em um único valor esperado corrigido — portanto, tratamos todo o fluxo de trabalho relacionado aos coeficientes e à combinação como pós-processamento.

Para avaliar até que ponto o MPF reflete a dinâmica real, calculamos primeiro os valores esperados exatos ao longo do tempo, elevando diretamente o hamiltoniano à potência de um. Isso só é viável porque $L = 10$; no exemplo de hardware em grande escala a seguir, teremos que recorrer a estimativas de redes tensoriais.

exact_expvals = []
for t in exact_evolution_times:
    exp_H = expm(-1j * t * hamiltonian.to_matrix())
    initial_state = Statevector(initial_state_circ).data
    time_evolved_state = exp_H @ initial_state

    exact_obs = (
        time_evolved_state.conj()
        @ observable.to_matrix()
        @ time_evolved_state
    ).real
    exact_expvals.append(exact_obs)

Coeficientes MPF estáticos

Os MPFs estáticos utilizam coeficientes $x_j$ que são independentes do tempo de evolução, do hamiltoniano e do estado inicial. Estabelecemos o sistema linear $Ax = b$ descrito na seção “Contexto” e calculamos os coeficientes. A matriz $A$ é determinada pelo número de passos de Trotter $k_j$, pela ordem $\chi$ da fórmula do produto e pelo fato de a fórmula ser simétrica (o que determina os expoentes $\eta_n$ ).

Para o nosso exemplo em pequena escala, utilizamos $k_j = [1, 2, 4]$ com uma fórmula de Suzuki-Trotter de ordem não simétrica — $2\chi=2$ (portanto, $\chi=1$ e $\eta_n = 2 + n$, o que resulta em $\eta_0 = 2,\, \eta_1 = 3$ ). O sistema passa a ser:

A primeira linha garante a imparcialidade ( $\sum_j x_j = 1$ ); a segunda e a terceira linhas cancelam, respectivamente, os termos de erro de Trotter de ordem principal $1/k^2$ e de ordem seguinte $1/k^3$.

Configure o LSE

Utilizamos setup_static_lse de qiskit_addon_mpf.static para construir a matriz $A$ e o vetor do lado direito $b$ descritos acima. A matriz $A$ depende não apenas de $k_j$, mas também da nossa escolha da fórmula do produto — em particular, sua ordem $\chi$ e se ela é simétrica. O symmetric sinalizador controla o padrão do expoente $\eta_n$ (fórmulas simétricas produzem apenas termos de erro de Trotter de potência par; ver Ref. [1] ). Observe que, conforme mostrado na Ref. [2], definir symmetric=True não é estritamente necessário, mesmo quando a função de potencial subjacente é simétrica — a LSE não simétrica continua válida (ela impõe restrições adicionais desnecessárias).

Para o nosso exemplo, já definimos order = 2 e symmetric = False na Etapa 1.

lse = setup_static_lse(mpf_trotter_steps, order=order, symmetric=symmetric)

Verifique a matriz construída $A$ e o vetor $b$ para confirmar se eles correspondem ao sistema descrito acima.

lse.A

Output:

array([[1.      , 1.      , 1.      ],
       [1.      , 0.25    , 0.0625  ],
       [1.      , 0.125   , 0.015625]])

lse.b

Output:

array([1., 0., 0.])

Com a equação de LSE em mãos, calculamos os coeficientes estáticos $x_j$ por meio de lse.solve() (essa é a solução direta $x = A^{-1}b$ ).

mpf_coeffs = lse.solve()
print(
    f"The static coefficients associated with the ansatze are: {mpf_coeffs}"
)

Output:

The static coefficients associated with the ansatze are: [ 0.04761905 -0.57142857  1.52380952]

Otimize para $x$ usando um modelo exato

Como alternativa ao cálculo de $x = A^{-1}b$, você pode usar setup\_exact\_model para construir uma instância de cvxpy.Problem que utilize o LSE como restrições e cuja solução ótima resulte em $x$.

model_exact, coeffs_exact = setup_exact_problem(lse)
model_exact.solve()
print(coeffs_exact.value)

Output:

[ 0.04761905 -0.57142857  1.52380952]

print(
    "L1 norm of the exact coefficients:",
    np.linalg.norm(coeffs_exact.value, ord=1),
)

Output:

L1 norm of the exact coefficients: 2.1428571428556378

Otimize para $x$ usando um modelo aproximado

Pode acontecer que a norma “ $L_1$ ” para o conjunto escolhido de valores de “ $k_j$ ” seja considerada muito alta. Se for esse o caso e você não puder escolher um conjunto diferente de valores de $k_j$, é possível usar uma solução aproximada que restrinja a norma $L_1$ a um limite escolhido, ao mesmo tempo em que minimiza $\|Ax - b\|$. Confira o guia sobre Como usar o modelo aproximado.

model_approx, coeffs_approx = setup_sum_of_squares_problem(
    lse, max_l1_norm=1.5
)
model_approx.solve()
print(coeffs_approx.value)
print(
    "L1 norm of the approximate coefficients:",
    np.linalg.norm(coeffs_approx.value, ord=1),
)

Output:

[-1.10294118e-03 -2.48897059e-01  1.25000000e+00]
L1 norm of the approximate coefficients: 1.5

Coeficientes dinâmicos do MPF

O MPF estático cancela os termos de erro de Trotter de uma forma independente do hamiltoniano e do estado; portanto, ele não produz necessariamente o menor erro de aproximação possível para um determinado hamiltoniano e estado inicial. O MPF dinâmico (Ref. [2], [3] ), por outro lado, determina coeficientes dependentes do tempo $x_i(t)$ que minimizam a distância na norma de Frobenius $\|\rho(t) - \mu^D(t)\|_F^2$ em cada instante $t$. Conforme mostrado na seção “Contexto”, isso requer a matriz de sobreposição $M_{ij}(t)$ entre os estados evoluídos por Trotter e a sobreposição $L_i(t)$ com o estado exato — ambos estimados por meio de backends de rede tensorial ( TeNPy ) em qiskit_addon_mpf.

Para configurar o LSE dinâmico, precisamos de três elementos:

1. Uma fábrica de evolutores aproximada que o complemento executará para cada $k_j$ a fim de produzir $\rho_{k_j}(t)$ como um MPS/MPO. Nós o construímos a partir da estrutura em camadas do circuito de Trotter de ordem $2$ (uma camada por slice_by_depth), envolvido como um LayerwiseEvolver com parâmetros de truncamento de TeNPy.
2. Uma fábrica de evolutores exatos que produz um $\rho(t)$ de referência de alta precisão. Utilizamos um circuito de Suzuki-Trotter de quarta ordem com pequeno passo de tempo (dt=0.1, order=4) como substituto da evolução exata.
3. Uma fábrica de identidades e um MPS de estado inicial que inicializam a simulação do TeNPy.

A célula abaixo cria a fábrica de evolutores aproximados.

# Create approximate time-evolution circuits
single_2nd_order_circ = generate_time_evolution_circuit(
    hamiltonian, time=1.0, synthesis=SuzukiTrotter(reps=1, order=order)
)
single_2nd_order_circ = pm.run(single_2nd_order_circ)  # collect XX and YY

# Find layers in the circuit
layers = slice_by_depth(single_2nd_order_circ, max_slice_depth=1)

# Create tensor network models
models = [
    LayerModel.from_quantum_circuit(layer, conserve="Sz") for layer in layers
]

# Create the time-evolution object
approx_factory = partial(
    LayerwiseEvolver,
    layers=models,
    options={
        "preserve_norm": False,
        "trunc_params": {
            "chi_max": 64,
            "svd_min": 1e-8,
            "trunc_cut": None,
        },
        "max_delta_t": 2,
    },
)

dt=0.1Aproximamos o estado exato ao longo do tempo com uma fórmula de Suzuki-Trotter de quarta ordem, utilizando um pequeno intervalo de tempo. Os parâmetros de truncamento do TeNPy podem afetar a precisão; portanto, é importante explorar uma variedade de valores.

single_4th_order_circ = generate_time_evolution_circuit(
    hamiltonian, time=1.0, synthesis=SuzukiTrotter(reps=1, order=4)
)
single_4th_order_circ = pm.run(single_4th_order_circ)
exact_model_layers = [
    LayerModel.from_quantum_circuit(layer, conserve="Sz")
    for layer in slice_by_depth(single_4th_order_circ, max_slice_depth=1)
]

exact_factory = partial(
    LayerwiseEvolver,
    layers=exact_model_layers,
    dt=0.1,
    options={
        "preserve_norm": False,
        "trunc_params": {
            "chi_max": 64,
            "svd_min": 1e-8,
            "trunc_cut": None,
        },
        "max_delta_t": 2,
    },
)

Por fim, definimos um identity_factory que gera o estado inicial do MPO e preparamos o estado inicial de Néel como um MPS que corresponde à rede utilizada pelo modelo de Trotter em camadas.

def identity_factory():
    return MPOState.initialize_from_lattice(models[0].lat, conserve=True)


mps_initial_state = MPS_neel_state(models[0].lat)

Com as fábricas definidas, calculamos agora os coeficientes dinâmicos em cada momento de evolução. Para cada $t$, setup_dynamic_lse calcula as matrizes de sobreposição relevantes por meio de TeNPy, e setup_frobenius_problem retorna um cvxpy.Problem que minimiza o custo da norma de Frobenius. O solucionador retorna coeficientes $x_j(t)$ específicos para esse momento; nós os reunimos em mpf_dynamic_coeffs_list. Se o solucionador falhar para um determinado $t$, voltamos a usar coeficientes iguais a zero para que o loop continue.

mpf_dynamic_coeffs_list = []
for t in trotter_times:
    print(f"Computing dynamic coefficients for time={t}")
    lse = setup_dynamic_lse(
        mpf_trotter_steps,
        t,
        identity_factory,
        exact_factory,
        approx_factory,
        mps_initial_state,
    )
    problem, coeffs = setup_frobenius_problem(lse)
    try:
        problem.solve()
        mpf_dynamic_coeffs_list.append(coeffs.value)
    except Exception as error:
        mpf_dynamic_coeffs_list.append(np.zeros(len(mpf_trotter_steps)))
        print(error, "Calculation Failed for time", t)
    print("")

Output:

Computing dynamic coefficients for time=0.5

Computing dynamic coefficients for time=0.6

Computing dynamic coefficients for time=0.7

Computing dynamic coefficients for time=0.7999999999999999

Computing dynamic coefficients for time=0.8999999999999999

Computing dynamic coefficients for time=0.9999999999999999

Computing dynamic coefficients for time=1.0999999999999999

Computing dynamic coefficients for time=1.1999999999999997

Computing dynamic coefficients for time=1.2999999999999998

Computing dynamic coefficients for time=1.4

Computing dynamic coefficients for time=1.4999999999999998

Combinar os valores esperados de Trotter com os coeficientes do MPF

Agora, calculamos o valor de “ $\langle A \rangle_{\text{MPF}}(t) = \sum_j x_j \, \langle A \rangle_{k_j}(t)$ ” para cada conjunto de coeficientes (estático-exato, estático-aproximado e dinâmico), propagamos os erros-padrão por circuito e representamos graficamente as séries temporais resultantes em relação à curva de diagonalização exata.

sym = {1: "^", 2: "s", 4: "p"}
# Get expectation values at all times for each Trotter step
for k, step in enumerate(mpf_trotter_steps):
    trotter_curve, trotter_curve_error = [], []
    for trotter_expvals, trotter_stds in zip(
        mpf_expvals_all_times, mpf_stds_all_times
    ):
        trotter_curve.append(trotter_expvals[k])
        trotter_curve_error.append(trotter_stds[k])

    plt.errorbar(
        trotter_times,
        trotter_curve,
        yerr=trotter_curve_error,
        alpha=0.5,
        markersize=4,
        marker=sym[step],
        color="grey",
        label=f"{mpf_trotter_steps[k]} Trotter steps",
    )

# Get expectation values at all times for the static MPF with exact coeffs
exact_mpf_curve, exact_mpf_curve_error = [], []
for trotter_expvals, trotter_stds in zip(
    mpf_expvals_all_times, mpf_stds_all_times
):
    mpf_std = np.sqrt(
        sum(
            [
                (coeff**2) * (std**2)
                for coeff, std in zip(coeffs_exact.value, trotter_stds)
            ]
        )
    )
    exact_mpf_curve_error.append(mpf_std)
    exact_mpf_curve.append(trotter_expvals @ coeffs_exact.value)

plt.errorbar(
    trotter_times,
    exact_mpf_curve,
    yerr=exact_mpf_curve_error,
    markersize=4,
    marker="o",
    label="Static MPF - Exact",
    color="purple",
)


# Get expectation values at all times for the static MPF with approximate coeffs
approx_mpf_curve, approx_mpf_curve_error = [], []
for trotter_expvals, trotter_stds in zip(
    mpf_expvals_all_times, mpf_stds_all_times
):
    mpf_std = np.sqrt(
        sum(
            [
                (coeff**2) * (std**2)
                for coeff, std in zip(coeffs_approx.value, trotter_stds)
            ]
        )
    )
    approx_mpf_curve_error.append(mpf_std)
    approx_mpf_curve.append(trotter_expvals @ coeffs_approx.value)

plt.errorbar(
    trotter_times,
    approx_mpf_curve,
    yerr=approx_mpf_curve_error,
    markersize=4,
    marker="o",
    label="Static MPF - Approx",
    color="orange",
)


# Get expectation values at all times for the dynamic MPF
dynamic_mpf_curve, dynamic_mpf_curve_error = [], []
for trotter_expvals, trotter_stds, dynamic_coeffs in zip(
    mpf_expvals_all_times, mpf_stds_all_times, mpf_dynamic_coeffs_list
):
    mpf_std = np.sqrt(
        sum(
            [
                (coeff**2) * (std**2)
                for coeff, std in zip(dynamic_coeffs, trotter_stds)
            ]
        )
    )
    dynamic_mpf_curve_error.append(mpf_std)
    dynamic_mpf_curve.append(trotter_expvals @ dynamic_coeffs)

plt.errorbar(
    trotter_times,
    dynamic_mpf_curve,
    yerr=dynamic_mpf_curve_error,
    markersize=4,
    marker="o",
    label="Dynamic MPF",
    color="pink",
)


# Exact expectation values
plt.plot(
    exact_evolution_times,
    exact_expvals,
    color="red",
    linestyle="--",
    label="Exact time-evolution",
)

plt.title(f"$\\langle Z_{{{L//2-1}}} Z_{{{L//2}}} \\rangle$ vs time")
plt.xlabel("Time")
plt.ylabel("Expectation Value")
plt.legend(loc="upper center", bbox_to_anchor=(0.5, -0.2), ncol=2)
plt.grid(alpha=0.1)
plt.tight_layout()
plt.show()

Output:

O gráfico acima ilustra a interação entre o erro de Trotter e o erro amostral.

• Erro de Trotter. As fórmulas individuais dos produtos (marcadores cinza) se afastam cada vez mais da curva exata à medida que o tempo passa. O circuito $k=1$ apresenta o maior desvio e é o mais raso, mas também já se encontra no regime em que $t/k \gtrsim 1$, de modo que o principal termo de erro $1/k^{2}$ é grande. As combinações de MPF (marcadores coloridos) anulam vários desses termos de erro de Trotter principais, de modo que acompanham a curva exata com muito mais precisão do que qualquer circuito de “ $k_j$ ” isolado. A lacuna restante reflete os termos de Trotter de ordem superior que o MPF não cancela: um MPF estático de ordem $2$, $r=3$ elimina apenas as duas primeiras ordens de erro e, em $t/k_{\min}$, a cauda não cancelada acaba por se tornar dominante — portanto, o MPF não garante que circuitos muito rasos permaneçam precisos em momentos arbitrários.

• Erro amostral. As barras de erro mais largas nas curvas do MPF são uma consequência direta da combinação linear: a propagação dos erros-padrão independentes por circuito $\sigma_{k_j}$ resulta em uma variância total $\sigma_{\text{MPF}}^2 = \sum_j x_j^2 \, \sigma_{k_j}^2$. Portanto, quanto maior for o $\|x\|_2$ (e, na prática, o $\|x\|_1$, que é o que controlamos), mais medições serão necessárias para atingir uma determinada incerteza alvo. Essa é a compensação por trás da opção “solucionador aproximado” na seção “Contexto”: limitamos o valor de $\|x\|_1$ para manter essa sobrecarga dentro de limites razoáveis. É fundamental ressaltar que, ao contrário do erro de Trotter, o erro amostral diminui com o aumento d $1/\sqrt{N_{\text{shots}}}$; portanto, ele sempre pode ser reduzido com o aumento do número de ensaios.

No exemplo de hardware em grande escala abaixo, o ruído do hardware surge como uma fonte adicional de erro em cada $\langle A \rangle_{k_j}$, que, da mesma forma, é amplificado pelos coeficientes do MPF. Veremos, nessa seção, como a mitigação de erros interage com os MPFs.

Exemplo de hardware em grande escala

Nesta seção, ampliemos a escala do problema para além do que é possível simular com exatidão. Reproduzimos alguns dos resultados apresentados na Ref. [3], utilizando uma cadeia XXZ de 50 qubits no tempo $t = 3$. Seguimos o mesmo fluxo de trabalho de quatro etapas do exemplo em pequena escala, agora voltado para hardware quântico real com mitigação de erros. Assim como no modelo, cada etapa é marcada diretamente no código, e uma única etapa pode abranger várias células quando vale a pena examinar seus resultados intermediários.

O mapeamento reflete o exemplo em pequena escala: definir um hamiltoniano, escolher os parâmetros de Trotter, calcular os coeficientes MPF (estáticos e dinâmicos) e construir circuitos. As principais diferenças são:

• Um hamiltoniano XXZ em 50 nós com acoplamentos aleatórios extraídos de $\mathcal{U}(0.5, 1.5)$ (Ref. [3] ).
• Uma fórmula de Trotter simétrica de segunda ordem com $k_j = [3, 4, 6]$ o (portanto, $\chi=1$, symmetric=True).
• Um único tempo de evolução fixo $t = 3$. Com $k_{\min}=3$, isso resulta em $t/k_{\min}=1$, mantendo os constituintes de profundidade rasa dentro do regime de convergência de Trotter, onde o modelo de erro dominante no qual o MPF se baseia é válido.
• Uma comparação adicional de circuito único foi executada com etapas de Trotter $k = 10$, usada como linha de base. Escolhemos o $k = 10$ porque sua profundidade de dois qubits no hardware é maior do que a do constituinte MPF mais profundo ( $k_{\max}=6$ ) somada à sobrecarga de executar múltiplos circuitos MPF — profundidade suficiente para ser limitada pelo ruído, que é o regime no qual se espera que a combinação de MPFs supere o desempenho da linha de base de circuito único. Trata-se de uma comparação de “circuito único profundo” com a combinação MPF, e não de um circuito que tenha como alvo o erro de Trotter efetivo do MPF (o que exigiria muito mais etapas).

Observe que, embora ainda estejamos na Etapa 1 (mapeamento e construção do circuito), também pré-calculamos os coeficientes dinâmicos juntamente com os estáticos nesta célula. Os coeficientes dinâmicos dependem de $H$ e $t$, mas não das medições quânticas; portanto, podem ser calculados a qualquer momento antes da Etapa 4. Fazemos isso agora para manter todas as configurações específicas do MPF em um único lugar.

# -------------------------Step 1-------------------------
L = 50
coupling_map = CouplingMap.from_line(L, bidirectional=False)

# XXZ Hamiltonian with random couplings (Ref. [3])
np.random.seed(0)
even_edges = list(coupling_map.get_edges())[::2]
odd_edges = list(coupling_map.get_edges())[1::2]

Js = np.random.uniform(0.5, 1.5, size=L)
hamiltonian = SparsePauliOp(Pauli("I" * L))
for i, edge in enumerate(even_edges + odd_edges):
    hamiltonian += SparsePauliOp.from_sparse_list(
        [
            ("XX", (edge), 2 * Js[i]),
            ("YY", (edge), 2 * Js[i]),
            ("ZZ", (edge), 4 * Js[i]),
        ],
        num_qubits=L,
    )

observable = SparsePauliOp.from_sparse_list(
    [("ZZ", (L // 2 - 1, L // 2), 1.0)], num_qubits=L
)

total_time = 3
mpf_trotter_steps = [3, 4, 6]
order = 2
symmetric = True

# Static coefficients
lse = setup_static_lse(mpf_trotter_steps, order=order, symmetric=symmetric)
mpf_coeffs = lse.solve()
print(f"Static coefficients: {mpf_coeffs}")
print(f"L1 norm: {np.linalg.norm(mpf_coeffs, ord=1)}")

model_approx, coeffs_approx = setup_sum_of_squares_problem(
    lse, max_l1_norm=2.0
)
model_approx.solve()
print(f"Approximate coefficients: {coeffs_approx.value}")
print(f"L1 norm (approx): {np.linalg.norm(coeffs_approx.value, ord=1)}")

# -------------------------Dynamic coefficients-------------------------
single_2nd_order_circ = generate_time_evolution_circuit(
    hamiltonian, time=1.0, synthesis=SuzukiTrotter(reps=1, order=order)
)
single_2nd_order_circ = pm.run(single_2nd_order_circ)

layers = slice_by_depth(single_2nd_order_circ, max_slice_depth=1)
models = [
    LayerModel.from_quantum_circuit(layer, conserve="Sz") for layer in layers
]

approx_factory = partial(
    LayerwiseEvolver,
    layers=models,
    options={
        "preserve_norm": False,
        "trunc_params": {"chi_max": 64, "svd_min": 1e-8, "trunc_cut": None},
        "max_delta_t": 4,
    },
)

single_4th_order_circ = generate_time_evolution_circuit(
    hamiltonian, time=1.0, synthesis=SuzukiTrotter(reps=1, order=4)
)
single_4th_order_circ = pm.run(single_4th_order_circ)
exact_model_layers = [
    LayerModel.from_quantum_circuit(layer, conserve="Sz")
    for layer in slice_by_depth(single_4th_order_circ, max_slice_depth=1)
]

exact_factory = partial(
    LayerwiseEvolver,
    layers=exact_model_layers,
    dt=0.1,
    options={
        "preserve_norm": False,
        "trunc_params": {"chi_max": 64, "svd_min": 1e-8, "trunc_cut": None},
        "max_delta_t": 3,
    },
)


def identity_factory():
    return MPOState.initialize_from_lattice(models[0].lat, conserve=True)


mps_initial_state = MPS_neel_state(models[0].lat)

print(f"Computing dynamic coefficients for time={total_time}")
lse_dyn = setup_dynamic_lse(
    mpf_trotter_steps,
    total_time,
    identity_factory,
    exact_factory,
    approx_factory,
    mps_initial_state,
)
problem, coeffs_dyn = setup_frobenius_problem(lse_dyn)
try:
    problem.solve()
    mpf_dynamic_coeffs = coeffs_dyn.value
except Exception as error:
    mpf_dynamic_coeffs = np.zeros(len(mpf_trotter_steps))
    print(error, "Calculation Failed")

# -------------------------Step 1 (cont): Build circuits-------------------------
mpf_circuits = []
for k in mpf_trotter_steps:
    circuit = QuantumCircuit(L)
    circuit.x([i for i in range(L) if i % 2])
    trotter_circ = generate_time_evolution_circuit(
        hamiltonian,
        synthesis=SuzukiTrotter(reps=k, order=order),
        time=total_time,
    )
    circuit.compose(trotter_circ, qubits=range(L), inplace=True)
    mpf_circuits.append(circuit)

# Baseline "single deep circuit" comparison run with k=10 Trotter steps.
# Its two-qubit depth is deeper than the deepest MPF constituent (k_max=6) plus
# the overhead of running multiple circuits, pushing it into the noise-limited
# regime where MPF is expected to outperform. It does NOT target the MPF's effective
# Trotter error (which would require many more steps).
comp_circuit = QuantumCircuit(L)
comp_circuit.x([i for i in range(L) if i % 2])
trotter_circ = generate_time_evolution_circuit(
    hamiltonian,
    synthesis=SuzukiTrotter(reps=10, order=order),
    time=total_time,
)
comp_circuit.compose(trotter_circ, qubits=range(L), inplace=True)
mpf_circuits.append(comp_circuit)

Output:

Static coefficients: [ 0.42857143 -1.82857143  2.4       ]
L1 norm: 4.65714285714286
Approximate coefficients: [-0.4942491   0.40206845  1.09218065]
L1 norm (approx): 1.9884981979026675
Computing dynamic coefficients for time=3

Agora otimizamos os circuitos para o backend escolhido. optimization_level=3Utilizamos o gerenciador de passagens predefinido do Qiskit, que seleciona automaticamente um bom conjunto de qubits físicos e direciona cada circuito para a topologia do dispositivo.

# -------------------------Step 2-------------------------
service = QiskitRuntimeService()
# backend = service.least_busy(operational=True, simulator=False, min_num_qubits=L)
backend = service.backend("ibm_fez")
print(backend)

transpiler = generate_preset_pass_manager(
    optimization_level=3, backend=backend
)
transpiled_circuits = [transpiler.run(circ) for circ in mpf_circuits]

isa_observables = [
    observable.apply_layout(circ.layout) for circ in transpiled_circuits
]

Output:

<IBMBackend('ibm_fez')>

A execução de circuitos mais complexos em hardware real exige medidas agressivas de mitigação de erros. Oferecemos descoplamento dinâmico, rotação de portas e medições, mitigação de erros de medição e extrapolação sem ruído (ZNE). Observe que os fatores de ruído ZNE que utilizamos aqui (1, 1.2, 1.4) são menores do que em um cenário de circuito raso, uma vez que os constituintes MPF mais profundos já estão próximos do limiar de ruído e grandes amplificações de ruído os levariam além do ponto em que a extrapolação ZNE é confiável.

Enviamos todos os quatro circuitos (três componentes do MPF em $k_j = [3, 4, 6]$, além da linha de base do $k = 10$ ) em uma única tarefa do Estimator.

# -------------------------Step 3-------------------------
estimator = Estimator(mode=backend)
estimator.options.default_shots = 30000

# Error suppression/mitigation
estimator.options.dynamical_decoupling.enable = True
estimator.options.twirling.enable_gates = True
estimator.options.twirling.enable_measure = True
estimator.options.twirling.num_randomizations = "auto"
estimator.options.twirling.strategy = "active-accum"
estimator.options.resilience.measure_mitigation = True
estimator.options.experimental.execution_path = "gen3-turbo"

estimator.options.resilience.zne_mitigation = True
estimator.options.resilience.zne.noise_factors = (1, 1.2, 1.4)
estimator.options.resilience.zne.extrapolator = "linear"

estimator.options.environment.job_tags = ["TUT_MPF"]

job_50 = estimator.run(
    [
        (circ, observable)
        for circ, observable in zip(transpiled_circuits, isa_observables)
    ]
)

Extraímos os valores esperados e os desvios-padrão por circuito dos resultados do trabalho e, em seguida, combinamo-los com cada conjunto de coeficientes MPF exatamente como no exemplo em pequena escala: $\langle A \rangle_{\text{MPF}} = \sum_j x_j \, \langle A \rangle_{k_j}$, com variância propagada $\sigma^2 = \sum_j x_j^2 \sigma_{k_j}^2$.

# -------------------------Step 4-------------------------
result = job_50.result()
evs = [res.data.evs for res in result]
std = [res.data.stds for res in result]

print(evs)
print(std)

Output:

[array(-0.07916195), array(-0.04479681), array(-0.2560756), array(-0.06045848)]
[array(0.04605538), array(0.10056336), array(0.14426151), array(0.04059092)]

exact_mpf_std = np.sqrt(
    sum([(coeff**2) * (std**2) for coeff, std in zip(mpf_coeffs, std[:3])])
)
print(
    "Exact static MPF expectation value: ",
    evs[:3] @ mpf_coeffs,
    "+-",
    exact_mpf_std,
)
approx_mpf_std = np.sqrt(
    sum(
        [
            (coeff**2) * (std**2)
            for coeff, std in zip(coeffs_approx.value, std[:3])
        ]
    )
)
print(
    "Approximate static MPF expectation value: ",
    evs[:3] @ coeffs_approx.value,
    "+-",
    approx_mpf_std,
)
dynamic_mpf_std = np.sqrt(
    sum(
        [
            (coeff**2) * (std**2)
            for coeff, std in zip(mpf_dynamic_coeffs, std[:3])
        ]
    )
)
print(
    "Dynamic MPF expectation value: ",
    evs[:3] @ mpf_dynamic_coeffs,
    "+-",
    dynamic_mpf_std,
)

Output:

Exact static MPF expectation value:  -0.5665938395816946 +- 0.3925273058119915
Approximate static MPF expectation value:  -0.25856647611537903 +- 0.164249927266166
Dynamic MPF expectation value:  -0.12667812062949296 +- 0.06059471006973169

sym = {3: "^", 4: "s", 6: "p"}
for k, step in enumerate(mpf_trotter_steps):
    plt.errorbar(
        k,
        evs[k],
        yerr=std[k],
        alpha=0.5,
        markersize=4,
        marker=sym[step],
        color="grey",
        label=f"{mpf_trotter_steps[k]} Trotter steps",
    )

plt.errorbar(
    3,
    evs[-1],
    yerr=std[-1],
    alpha=0.5,
    markersize=8,
    marker="x",
    color="blue",
    label="10 Trotter steps",
)

plt.errorbar(
    4,
    evs[:3] @ mpf_coeffs,
    yerr=exact_mpf_std,
    markersize=4,
    marker="o",
    color="purple",
    label="Static MPF",
)

plt.errorbar(
    5,
    evs[:3] @ coeffs_approx.value,
    yerr=approx_mpf_std,
    markersize=4,
    marker="o",
    color="orange",
    label="Approximate static MPF",
)

plt.errorbar(
    6,
    evs[:3] @ mpf_dynamic_coeffs,
    yerr=dynamic_mpf_std,
    markersize=4,
    marker="o",
    color="pink",
    label="Dynamic MPF",
)

exact_obs = -0.24384471447172074  # Calculated via Tensor Network calculation
plt.axhline(
    y=exact_obs, linestyle="--", color="red", label="Exact time-evolution"
)

plt.title(
    f"$\\langle Z_{{{L//2-1}}} Z_{{{L//2}}} \\rangle$ at time {total_time} for the different methods"
)
plt.xlabel("Method")
plt.ylabel("Expectation Value")
plt.legend(loc="upper center", bbox_to_anchor=(0.5, -0.2), ncol=2)
plt.grid(alpha=0.1)
plt.tight_layout()
plt.show()

Output:

Algumas observações sobre os resultados de hardware apresentados acima:

• Aprofundar-se não é gratuito em termos de hardware. As linhas de base de circuito único mostram claramente o que acontece: o circuito $k = 6$ é praticamente exato ( $-0.256$ em comparação com a referência $-0.244$ ), mas a linha de base mais profunda $k = 10$ apresenta um resultado pior ( $-0.061$, com um desvio de $\sim 0.18$ ), e não melhor. Quando o erro de Trotter já é pequeno, adicionar etapas basicamente aumenta a profundidade do circuito e acumula mais ruído de porta e decoerência. É exatamente para esse regime que os MPFs foram concebidos: alcançar a precisão de um circuito profundo utilizando apenas componentes superficiais.

• Um MPF de norma pequena supera o circuito único profundo. O MPF estático aproximado (limitado a $\|x\|_1 \approx 2$ ) fica em $-0.259$, a $\sim 0.015$ do valor de referência e muito mais próximo do que a linha de base de $k = 10$. O MPF dinâmico ( $-0.127$ ) também supera com folga essa referência. Ambos combinam apenas os circuitos rasos $k_j = [3, 4, 6]$, mas conseguem chegar a uma resposta que o circuito único profundo não conseguiu.

• A norma do coeficiente é mais importante do que a otimização matemática. O MPF estático exato possui um $\|x\|_1 = 4.66$ e e é o pior estimador de todos ( $-0.567$, com um desvio de mais de $0.3$ ): a grande norma do coeficiente amplifica o ruído residual do portão, a decoerência e o erro ZNE em cada $\langle A \rangle_{k_j}$ aproximadamente pelo mesmo fator, anulando o cancelamento do erro de Trotter que ele proporciona. A limitação da norma (o solucionador estático aproximado, $\|x\|_1 \approx 2$ ) elimina essa sobrecarga e fornece a melhor estimativa — mesmo que seus coeficientes não cancelem mais exatamente o erro principal de Trotter.

• Circuitos individuais de pouca profundidade ainda podem ser competitivos. O único componente do método de fluxo de massa ( $k = 6$ ) ( $-0.256$ ) é, por si só, essencialmente exato neste caso — nesta simulação, ele é até mesmo ligeiramente mais preciso do que o MPF estático aproximado. O problema é que você não sabe de antemão qual $k$ se encontra no ponto ideal de “convergência alcançada, mas ainda não limitada pelo ruído”, e a opção que parece segura — simplesmente ir mais fundo ( $k = 10$ ) para garantir a convergência de Trotter — é exatamente aquela que falha. O MPF oferece uma combinação baseada em princípios de circuitos rasos que não exige adivinhar a profundidade correta.

A lição prática é que, no hardware, os MPFs devem ser combinados com uma forte mitigação de erros em cada $\langle A \rangle_{k_j}$ individual; a norma do coeficiente $L_1$ deve ser mantida modesta (use o solucionador aproximado ou o MPF dinâmico); e os passos de Trotter $k_j$ devem ser escolhidos de modo que $t/k_{\min} \lesssim 1$ — aqui, $k_{\min} = 3$ em $t = 3$ resulta em $t/k_{\min} = 1$, mantendo os constituintes dentro do regime de convergência, onde o modelo de erro principal no qual o MPF estático se baseia é válido. Com essas escolhas, os MPFs de norma pequena aqui apresentados se equiparam a um circuito único convergente, ao passo que a linha de base ingênua do tipo “basta ir mais fundo” não o faz, recuperando a vantagem de profundidade versus precisão demonstrada na Ref. [3]. Observe também que as execuções individuais apresentam ruído — em um envio diferente do mesmo trabalho (ou em um backend diferente), a ordem exata pode mudar; as tendências consistentes são: os MPFs de pequeno porte $\|x\|_1$ apresentam bom desempenho, o MPF de grande porte $\|x\|_1$ com estática exata é amplificado pelo ruído do hardware e o circuito único excessivamente profundo é limitado pelo ruído.

Próximas etapas

Referências

[1] Vázquez, A. C., Egger, D. J., Ochsner, D., e Woerner, S. Fórmulas multiproduto bem condicionadas para simulação hamiltoniana otimizada para hardware. Quantum, 7, 1067 (2023)

[2] Zhuk, S., Robertson, N. F., & Bravyi, S. Limites de erro de Trotter e fórmulas dinâmicas multiproduto para simulação hamiltoniana. Physical Review Research, 6(3), 033309 (2024)

[3] Robertson, N. F., et al. Fórmulas dinâmicas multiproduto aprimoradas por redes tensoriais. arXiv:2407.17405 (2024)