I am an ITensor user who recently asked about custom Majorana fermions. I have a question regarding the calculation of four-point correlation functions for a system of spinless fermions using Matrix Product States (MPS).

I understand how to compute two-point correlation functions using the correlation_matrix function. However, I’m unsure how to proceed with four-point functions (for example, <c_1c_2c_3c_4>). In spin systems, I’ve learned to use the “apply” method from a senior in my lab, but I suspect this might not work the same way in fermionic systems due to the presence of Jordan-Wigner strings. Indeed, when I compared the results from the correlation_matrix method and the “apply” method for two-point functions, they yielded different results.

Could anyone advise on the correct approach to calculate four-point functions for fermionic systems using MPS in ITensor? Any guidance or reference to relevant documentation would be greatly appreciated.

Good question, so this kind of touches on a fundamental challenge for MPS and many-body methods. For example, people working with Green’s functions and Feynman diagrams also are working on new methods to handle “vertex functions” which involve four variables.

The reason it’s challenging is this: say we had a function which could instantly calculate C_{i_1 i_2 i_3 i_4} = \langle A_{i_1} B_{i_2} C_{i_3} D_{i_4} \rangle . If the number of lattice sites is N, just storing this would cost N^4 memory. For just N=100, this is a hundred million numbers (N^4 = 10^9). So you can imagine that analyzing this quantity is even harder, and it needs to be approached with care.

That being said, we do have an experimental package you can try: https://github.com/ITensor/ITensorCorrelators.jl
which computes four-point correlators in an efficient way, and should handle fermions automatically as long as the relevant operators are set up properly (i.e. come from inside the ITensor library and/or have the has_fermion_string function defined for them). I would start by trying this package on very small (N\sim10) systems and checking against exact results, e.g. free fermions or ED.

Another thing you can do is to think through your whole calculation pipeline from end to end. Do you really need the entire 4-point correlator C_{i_1 i_2 i_3 i_4} or just certain “slices” of it or certain linear combinations. Like maybe you are going to take specific Fourier transforms of it afterward. If so, then you can do tricks like write down the operators which are already the Fourier transform or, more generally, linear combination you want then make MPOs of these operators and just measure these MPOs. Like say you want to measure some “structure factors” or similar at specific Brillouin zone points (\Gamma point, K point, …). Then you could make MPOs that just measure each of these specific points only.

I understand the suggestion regarding the experimental package. In fact, what I actually need are specific “slices” of the four-point function (for example, <c_{i}c_{j}c_{j+1}c_{j+2}>), so I think I can manage the memory usage with some careful planning.

I will start by testing this approach in a smaller system as you advised. Thank you again for your help.

I’ve encountered an issue with the correlator function in ITensor when trying to calculate correlations on the same site, such as correlator(psi, ("Cdag", "Cdag"), [(500, 500)]). I am getting a DimensionMismatch error as follows:

vbnetCopy code

DimensionMismatch: In scalar(T) or T[], ITensor T is not a scalar (it has indices ((dim=101|id=747|"Link,l=499"), (dim=101|id=747|"Link,l=499")', (dim=101|id=455|"Link,l=500"), (dim=101|id=455|"Link,l=500")')).

Theoretically, I can describe what I want to calculate using the operator “N”, but it’s quite complex to implement. Is there a way to circumvent this error or a workaround that I might consider?