Accuracy of split-gate MPS state evolution

[vincentmr]

I slightly modified the example in the Quimb doc showing how to compute local expectation values in circuits with MPS states (see below). The doc says we have the following options when applying gates:

• contract=False
• contract='split-gate'
• contract=True
• contract='swap+split'

The way I understand it is the first two options are equivalent to building a tensor network describing the circuit. The third option contracts everything on-the-fly, so if one has a decent amount of entanglement, one quickly ends up with a 2**N tensor, pretty much like a wavefunction. The fourth and last one seems to correspond to evolving an MPS state. Like the doc says: "This explicitly maintains the exact structure of an MPS (at the cost of increasing bond dimension)."

When I increase entanglement in a system, all methods agree except swap+split. I would expect that if the bond dimension is free to grow, agreement with other methods should be maintained. Am I wrong?

# %% Import Quimb + defs
import quimb as qu
from quimb.tensor.circuit import MPS_computational_state
from quimb.tensor.tensor_builder import MPO_rand_herm
from quimb.gates import CNOT, RY, RZ, H

# %%
n = 10
A = MPO_rand_herm(n, bond_dim=2, tags=["HAM"])


# %%
def two_qubit_layer(circ, gate_round, contract=False):
    """Apply a layer of constant entangling gates."""
    n_wires = circ.L
    offset = (gate_round % (n_wires - 1)) + 1
    for i in range(circ.L):
        circ.gate_(CNOT, (i, (i + offset) % n_wires), contract=contract)


def circuit(n_layers, contract=False):
    psi0 = MPS_computational_state("0" * n)

    for j in range(n_layers):
        for i in range(n):
            psi0.gate_(H, i, contract=contract)

        Rz = qu.phase_gate(0.42)
        for i in range(n):
            psi0.gate_(Rz, i, contract=contract)

        two_qubit_layer(psi0, j, contract=contract)
    return psi0


# %%
for nl in range(1, 4):
    for contract in [False, "split-gate", True, "swap+split"]:
        circ = circuit(nl, contract=contract)

        p = circ
        pH = circ.H
        p.align(A, pH)
        expv_qub = (pH & A & p).contract(all, output_inds=())
        print(expv_qub)

# (1.687574551021768e-05+3.012594425061628e-22j)
# (1.68757455102175e-05+9.172034541602171e-21j)
# (1.6875745510217717e-05+0j)
# 1.6875745510217877e-05
# (0.0011250328588973567+6.776263578034403e-20j)
# (0.001125032858897338-5.909775884721955e-21j)
# (0.001125032858897356+0j)
# (-0.0014766939810730935+0j)
# (-0.00019349567437897374-1.6940658945086007e-20j)
# (-0.00019349567437896908-4.319868030996932e-20j)
# (-0.0001934956743789734+0j)
# (-0.0001398444705628218+0j)

[jcmgray]

The mps methods by default use a small but non-zero cutoff, do you get the same results if you set cutoff=0.0, max_bond=None which should be the exact equivalent?

[vincentmr]

Do you mean modifying the script as

gate_opts = dict(cutoff=0.0, max_bond=None)


# %%
def two_qubit_layer(circ, gate_round, contract=False):
    """Apply a layer of constant entangling gates."""
    n_wires = circ.L
    offset = (gate_round % (n_wires - 1)) + 1
    for i in range(circ.L):
        circ.gate_(CNOT, (i, (i + offset) % n_wires), contract=contract, **gate_opts)


def circuit(n_layers, contract=False):
    psi0 = MPS_computational_state("0" * n)

    for j in range(n_layers):
        for i in range(n):
            psi0.gate_(H, i, contract=contract, **gate_opts)

        Rz = qu.phase_gate(0.42)
        for i in range(n):
            psi0.gate_(Rz, i, contract=contract, **gate_opts)

        two_qubit_layer(psi0, j, contract=contract)
    return psi0

? It yields similar results. I appears cutoff=0.0, max_bond=None are not arguments of the MPS_computational_state constructor (or MatrixProductState), is that right?

[jcmgray]

I think the issue might actually be a bug where gate_with_auto_swap treats the gate being applied as symmetric and calls i, j = sorted(where). Could you try with gates where i < j always to check?

[jcmgray]

Have pushed a fix here - b91f7a6 - again, let me know if that works!

[vincentmr]

Right on, everything matches perfectly now! Thanks for your continued support. I'll mark your original answer as the answer, alright?