Mismatch between `correlator` and `expectation_value` for fermionic PEPS

[Yue-Zhengyuan]

It seems that the correlation function calculated from correlator (added in #210) is different from expectation_value for fermionic PEPS.

To reproduce the issue, run the following code that calculates the spin correlation $\langle S_i \cdot S_j \rangle$ along the horizontal direction in t-J model using a random iPEPS and a random CTMRGEnv (I'm now using TensorKitTensors instead of MPSKitModels):

using Random
using TensorKit
using PEPSKit
import TensorKitTensors: TJOperators as tJ

Pspace = tJ.tj_space(Trivial, Trivial)
Vspace = Vect[FermionParity](0 => 2, 1 => 2)
Espace = Vect[FermionParity](0 => 3, 1 => 3)
Random.seed!(100)
peps = InfinitePEPS(rand, Float64, Pspace, Vspace; unitcell=(2, 2));
env = CTMRGEnv(rand, Float64, peps, Espace);
lattice = collect(space(t, 1) for t in peps.A)

site0 = CartesianIndex(1, 1)
maxsep = 6
site1s = collect(site0 + CartesianIndex(0, i) for i in 1:maxsep)
op = tJ.S_exchange(Float64, Trivial, Trivial);

vals1 = collect(
    begin
        O = LocalOperator(lattice, (site0, site1) => op)
        val = expectation_value(peps, O, env)
    end for site1 in site1s
)
vals2 = correlator(peps, op, site0, site1s, env)
@show @show vals1 ≈ vals2 # false

But if I calculate the spin correlation in Heisenberg model, correlator and expectation_value produce (approximately) the same results. This is why the test test/examples/tf_ising.jldidn't catch the issue. I think we need to change how the column transfer matrix acts on the already-contracted part of the correlator "sandwich": for fermionic iPEPS, simply using * may have missed some twists.

PEPSKit.jl/src/algorithms/correlators.jl

Lines 39 to 41 in 47e6675

	T = TransferMatrix(Atop, sandwich, _dag(Abot))
	Vo = Vo * T
	Vn = Vn * T

[lkdvos]

Yeah, I also have noticed this problem at some point. The * is actually calling out to MPSKit.transfer_left, but unfortunately there is a missing twist in that thing. Most likely this is due to the interplay between the _dag to undo that this function puts an additional conj in the contraction, but we'd have to investigate a bit more to see what's going on

[Yue-Zhengyuan]

https://github.com/QuantumKitHub/MPSKit.jl/blob/2a7856c68eb9c1a79ee251b4439b68f120a7eb5d/src/transfermatrix/transfer.jl#L18-L27

I see that @plansor is used there. I guess it's because MPS calculations regards the transfer matrix as a linear operator. Here we can use twist! on Vo, Vn if we want to avoid modifying MPSKit.

[lkdvos]

They are overloaded though, using @tensor here:

https://github.com/QuantumKitHub/PEPSKit.jl/blob/master/src/algorithms/contractions/vumps_contractions.jl#L43

[Yue-Zhengyuan]

I have some time this week. I'll try if I can fix it.