Expectation values

Expectation values in AKLT state

def measure_mps_aklt(config=None, tol=1e-12):
    """ Test measuring MPS expectation values in AKLT state"""
    #
    # AKLT state with open boundary conditions and N = 31 sites.
    #
    N = 31
    psi = build_aklt_state_manually(N=N, lvec=(1, 0), rvec=(1, 0), config=config)
    #
    # We verify transfer matrix in the middle of AKLT state
    # to specify the effect of lvec and rvec
    #
    A = psi[N // 2]  # read MPS tensor
    T = yastn.ncon((A, A.conj()), ((-0, 1, -2), (-1, 1, -3)))
    T = T.fuse_legs(axes=((0, 1), (2, 3)))
    Tref = np.array([[1,  0,  0,  2],
                     [0, -1,  0,  0],
                     [0,  0, -1,  0],
                     [2,  0,  0,  1]])
    assert np.allclose(T.to_numpy(), Tref / 3, atol=tol)
    #
    # psi1 is not normalized.
    # Norm is very close to sqrt(0.5) for this system size.
    #
    norm = psi.norm()  # use canonization to calculate norm
    norm2 = mps.vdot(psi, psi)  # directly contract mps squared
    #
    assert abs(norm - np.sqrt(0.5)) < tol
    assert abs(norm2 - 0.5) < tol
    #
    # Initialize dense Spin1 operators for expectation values calculations
    #
    opts_config = {} if config is None else \
                  {'backend': config.backend,
                   'default_device': config.default_device}
    # pytest uses config to inject various backends and devices for testing
    ops = yastn.operators.Spin1(sym='dense', **opts_config)
    #
    # Measure Sz at the first and last sites
    #
    ez = mps.measure_1site(psi, {0: ops.sz(), N-1: ops.sz()}, psi)
    assert len(ez) == 2 and isinstance(ez, dict)
    assert abs(ez[0] - 1./3) < tol  # psi is not normalized
    assert abs(ez[N - 1] + 1./3) < tol
    #
    #  the same with other syntaxes
    #
    ez_bis = mps.measure_1site(psi, ops.sz(), psi, sites=[0, N-1])
    assert ez == ez_bis
    ez0 = mps.measure_1site(psi, ops.sz(), psi, sites=0)  # this return a number
    assert abs(ez0 - 1./3) < tol
    #
    #  with normalised state
    psin = psi / norm
    ez = mps.measure_1site(psin, ops.sz(), psin)  # calculate sz on all sites
    assert len(ez) == N and isinstance(ez, dict)
    assert abs(ez[0] - 2. / 3) < tol  # psin is normalized
    #
    # Calculate <Sz_i Sz_j> for (i, j) = (n, n+r) with r=1
    ezz = mps.measure_2site(psin, ops.sz(), ops.sz(), psin, bonds='r1')
    assert len(ezz) == N - 1 and isinstance(ezz, dict)
    assert all(abs(ezznn + 4. / 9) < tol for ezznn in ezz.values())
    #
    # here calculate a single bond, returning a number
    ezz12 = mps.measure_2site(psin, ops.sz(), ops.sz(), psin, bonds=(1, 2))
    assert abs(ezz12 + 4. / 9) < tol
    #
    #  Measure string operator; exp (i pi * Sz) = I - 2 * abs(Sz)
    #
    esz = ops.I() - 2 * abs(ops.sz())
    Hterm = mps.Hterm(amplitude=1,
                      positions=(10, 11, 12, 13, 14, 15),
                      operators=(ops.sz(), esz, esz, esz, esz, ops.sz()))
    I = mps.product_mpo(ops.I(), N)
    Ostring = mps.generate_mpo(I, [Hterm])  # MPO a string correlator
    Estring = mps.vdot(psin, Ostring, psin)
    assert abs(Estring + 4. / 9) < tol

Entropy and spectrum in GHZ-like state

def mps_spectrum_ghz(sym='dense', config=None, tol=1e-12):
    """ Test measuring Schmidt spectrum and entropy of mps. """
    N = 8
    # consider 3-dimensional local Hilbert space
    opts_config = {} if config is None else \
                    {'backend': config.backend,
                    'default_device': config.default_device}
    # pytest uses config to inject various backends and devices for testing
    ops = yastn.operators.Spin1(sym=sym, **opts_config)
    #
    # take 3 orthonormal product states
    #
    psi1 = mps.product_mps(ops.vec_z(val=+0), N)
    psi2 = mps.product_mps([ops.vec_z(val=+1), ops.vec_z(val=-1)], N)
    psi3 = mps.product_mps([ops.vec_z(val=-1), ops.vec_z(val=+1)], N)
    #
    # For Z3 and U1, their sum has a well-defined
    # global charge only for even N
    #
    psi = mps.add(psi1, psi2, psi3, amplitudes=(0.5, 0.25, -0.25))
    #
    # state psi is not normalized
    #
    assert abs(psi.norm() - np.sqrt(0.375)) < tol
    #
    # calculate Schmidt values; Schmidt values are properly normalized
    #
    schmidt_spectra = psi.get_Schmidt_values()
    assert len(schmidt_spectra) == N + 1 and isinstance(schmidt_spectra, list)
    assert all(abs(sv.norm() - 1) < tol for sv in schmidt_spectra)
    #
    bds = (1,) + (3,) * (N-1) + (1,)  # (1, 3, 3, 3, 3, 3, 3, 3, 1) for N = 8
    assert psi.get_bond_dimensions() == bds
    # schmidt values are stored as diagonal yastn.Tensors
    assert all(sv.isdiag and sv.get_shape() == (bd, bd)
               for sv, bd in zip(schmidt_spectra, bds))
    #
    # mps.get_Schmidt_values is used by mps.get_entropy
    #
    entropies = psi.get_entropy()
    assert len(entropies) == N + 1 and isinstance(entropies, list)
    #
    # by default calculate von Neuman entropy (base-2 log)
    #
    assert all(abs(ent - (np.log2(6) - 4/3)) < tol for ent in entropies[1:-1])
    assert abs(entropies[0]) < tol and abs(entropies[-1]) < tol
    #
    # Renyi entropy with alpha=2 (and base-2 log)
    #
    entropies2 = psi.get_entropy(alpha=2)
    assert all(abs(ent - 1.) < tol for ent in entropies2[1:-1])
    assert abs(entropies2[0]) < tol and abs(entropies2[-1]) < tol
    #
    # the state psi is not affected by calculations
    #
    assert abs(psi.norm() - np.sqrt(0.375)) < tol