Fourier Transform of Hamiltonian and MPO

DMRG and Numerical Methods
, ,
1

Hello everyone,

I am currently working on extracting critical data, such as central charge and scaling dimension, from a 1D critical quantum spin chain. To achieve this, I am following the approach outlined in equation 8 of the paper https://arxiv.org/pdf/1706.01436.pdf which involves constructing a Fourier transform of the Hamiltonian.

image
image816×166 5.05 KB

I am interested in exploring the possibility of creating an MPO (Matrix Product Operator) version for the Fourier-transformed Hamiltonian within the ITensor library using Julia. More specifically, I am aiming to implement the expression H_2=\sum exp(4 piij/N) H_j from equation 8.

If anyone has insights, suggestions, or code snippets related to this particular task, I would greatly appreciate your input. Any guidance on how to effectively construct this MPO version for the Fourier-transformed Hamiltonian would be extremely helpful.

Thank you in advance for your time and assistance!

2

Hello,

I may be able to help you. I assume you already have your H0 as an matrix product operator. Then what you want to create is H2, is that the case?

Your operation \exp(ij4\pi/N) act site by site individually, thus you can define an operator and then apply it directly to H0 or create an MPO which can then be applied to H0:

#operator creation
function op(::OpName"QFT", ::SiteType"Qubit"; j, N)
    a = exp(1im*j*4π/N)
    return [
        a 0
        0 a
    ]
end

Then to apply this opperator I suggest two possible solutions:

N = 10                               #length of your circle
cut = 1e-12
s = siteinds("Qubit", N)     #site creation
H0 = MPO(s)                    #create H0
H2 = H0
for i=1:N
    G = op("QFT", s[i]; j=i,N=N)
    H2 = apply(G, H2; cutoff=cut)
end

or you can also do

N = 10                               #length of your circle
cut = 1e-12
s = siteinds("Qubit", N)     #site creation
H0 = MPO(s)                    #create H0
QFT_MPO = MPO(s, "Id")
for i=1:N
    G = op("QFT", s[i]; j=i,N=N)
    QFT_MPO = apply(G, QFT_MPO; cutoff=cut)
end
H2 = apply(QFT_MPO, H0; cutoff=cut)

Hope this helps,
Aleix

3

Hi Alex,

Thank you for your prompt reply. Since I have limited experience with custom MPOs, I have some follow-up questions to better understand how this custom MPO interacts with local Hamiltonians.

You provided the following MPO construction and apply it to H_0:

function op(::OpName"QFT", ::SiteType"Qubit"; j, N)
a = exp( i* 4π*j / N) 
return [a 0; 0 a] 
end

Could you elaborate on how this MPO acts on the local Hamiltonian h_j? I am particularly interested in this because my Hamiltonian H includes cubic terms also.

For instance, when using closed boundary conditions, my Hamiltonian is defined as:

for j = 1:N-1
    ampo += -t1/2 ,"Sx",j,"Sx",j+1
    ampo += J/4 ,"Sz",j,"Sz",j+1
end
ampo += -t1/2 ,"Sx",N,"Sx",1
ampo += J/4 ,"Sz",N,"Sz",1

for j = 1:N-2
    ampo += t2 ,"Sx",j,"Sz",j+1,"Sx",j+2
end
ampo += t2 ,"Sx",N-1,"Sz",N,"Sx",1
ampo += t2 ,"Sx",N,"Sz",1,"Sx",2

H = MPO(ampo,sites)

Alternatively, with open boundary conditions, it is:

ampo = OpSum()
for j = 1:N-1
    ampo += -t1/2 ,"Sx",j,"Sx",j+1
    ampo += J/4 ,"Sz",j,"Sz",j+1
end
for j = 1:N-2
    ampo += t2 ,"Sx",j,"Sz",j+1,"Sx",j+2
end
H = MPO(ampo,sites)

I would really appreciate any explanation and further insights you can offer on this.

Thank you so much!

4

Hello,

The op() that I have defined acts on a single site of your Hamiltonian H. You should look at the op() as a single qubit gate. You can always create an MPO of rank 1 as a product of one qubit gates, that is what I am doing in the second example. In the first example, I apply the one qubit gate directly on your Hamiltonian (the results should be the same).

From what I understand the Fourier Transform always have periodic boundary conditions. However, one key point that I would like to emphasize is that the Fourier Transform operator is independent of the Hamiltonian given. In this sense, it is similar to the Fourier Transform in the continuous space where
to transform a function into the fourier space you have to apply \int_x e^{iqx}.

I cannot comment how the different boundary conditions will affect your results. My guess is that probably it is better to work with a periodic (closed system) but I am unsure about that.

5

Hi @ajitkumar1991, if it’s helpful for your case, the OpSum system can accept arbitrary long-ranged terms and terms with arbitrary numerical prefactors.

So if e.g. you wanted to make an MPO for an operator like

O = \sum_{j=1}^{N-1} e^{i j} S^z_j S^z_{j+1}

then you could do this by writing:

os = OpSum()
for j=1:N-1
  os += exp(im*j), "Sz",j, "Sz", j+1
end
6

Hi Aleix,

Sorry for the delay in my response. I reviewed your suggestions and it is working as expected now.
Thank you so much!

1 Like
7

Hi Miles,

Thank you for your suggestion. I will try it and let you know if I have more questions.