Combiner built-in function issue for virtual bond index from MPO with quantum number conservation


Now, I am coding 1-site DMRG with subspace expansion method (or 3S method) in Julia, which includes the additional step of adding perturbation term by increasing virtual bond dimension. Initially, I encountered an issue during testing while using the combiner function to merge Indexes considering quantum numbers.
When selecting two virtual bond index each one from the MPS and MPO of a specific site and combining them with the combiner function, it appears that the dimensions of the quantum numbers do not match. Selecting two virtual bond indexes from the MPS and combining them with the combiner function results in properly matching quantum numbers and dimensions.
However, upon examining the indexes of the MPO, All virtual bond indices include empty tagged quantum number but assigned dimension 1 (QN() => 1), indicating that quantum numbers have not been assigned. Furthermore, when combining the indexes of the MPS and MPO, while some of quantum numbers correctly match the bond dimension, others yield larger dimensions. It’s unclear what the issue is.

(P.S. We want to perform a direct sum of tensors, but currently, we’re copying each component one by one using a for loop, which is too slow. Is there a different method available?)"

Could you share some minimal runnable example code demonstrating the issue you are seeing?

There is a directsum function available in ITensors.jl for direct summing two (or more) ITensors.

Sorry for late due to lunch

Test code is as follows.


sites = siteinds(“S=1”, N; conserve_qns= true )

os = OpSum()

for j=1:N-1

os += “Sz”,j,“Sz”,j+1

os += 1/2,“S+”,j,“S-”,j+1

os += 1/2,“S-”,j,“S+”,j+1


H = MPO(os, sites)

states = [isodd(n) ? “Up” : “Dn” for n in 1:N]

χ = [ i for i in 1:N-1]
psi = randomMPS(sites, states; linkdims=χ)
id_mps = inds(psi[6], tags=“Link,l=6”)[1]
id_mpo = inds(H[6], tags= “Link,l=6”)[1]
id_comb = combiner(id_mps, id_mpo; tags=“id_comb”)

This gives resulting that:
ITensor ord=3
1: QN(“Sz”,-4) => 2
2: QN(“Sz”,-2) => 9
3: QN(“Sz”,0) => 12
4: QN(“Sz”,2) => 6
5: QN(“Sz”,4) => 1
1: QN(“Sz”,2) => 1
2: QN(“Sz”,0) => 3
3: QN(“Sz”,-2) => 2
1: QN() => 1
2: QN(“Sz”,0) => 1
3: QN(“Sz”,2) => 1
4: QN(“Sz”,-2) => 1
5: QN(“Sz”,0) => 1

where question is three.

Q1. For id=152 index, why QN() => 1 exists?

Q2. For id=152 index, why the two quantum numbers are separated 2: QN(“Sz”,0) => 1 and 5: QN(“Sz”,0) => 1

Q3. Why the combinations of incoming indices (id=152, id=176) are not match to the resulting bond dimension of outgoing index (id=785)?
The example combinations are follows.

QN(Sz=-2, id=152) + QN(Sz=0, id=176) gives dimension => 2 * 1
QN(Sz=-2, id=152) + QN(Sz=0, id=176) gives dimension => 2 * 1
QN(Sz=0, id=152) + QN(Sz=-2, id=176) gives dimension => 3 * 1
It would be that the total bond dimension is 7 for QN(Sz=-2) outgoing.
But resulting bond dimension QN(Sz=-2, id=785) => 9.

When i command that directsum(A, B), it gives error: MethodError: no method matching directsum(::ITensor, ::ITensor)

(where A, B are 3-legs tensors that two indices are shared and a index are different)

Please look at the documentation for directsum by typing ?directum at the Julia prompt for instructions on how to use it, you have to specify which indices should get direct summed.

Please format your code, see the post Please Read: Make It Easier to Help You.

Q1. and Q2. That index is constructed as part of the OpSum to MPO conversion code. The code creates indices with that structure to keep track of operators in different QN sectors. It isn’t a problem that the sectors are split in that way, in fact it can lead to MPO tensors that are more sparse which leads to better DMRG performance. QN() is a “universal zero” sector so it will be treated as QN("Sz", 0) in this case.

Q3. The part you are missing is that the QN() => 1 is a “universal zero” sector, so is treated as another QN("Sz", 0) sector in this case.

Here is the running code

using ITensors

sites = siteinds("S=1", N; conserve_qns= true );
os = OpSum();

for j=1:N-1
    os += "Sz",j,"Sz",j+1
    os += 1/2,"S+",j,"S-",j+1
    os += 1/2,"S-",j,"S+",j+1

H = MPO(os, sites);

states = [isodd(n) ? "Up" : "Dn" for n in 1:N];
χ = [ i for i in 1:N-1];
psi = randomMPS(sites, states; linkdims=χ);
id_mps = (inds(psi[6], tags="Link,l=6")[1]);
id_mpo = (inds(H[6], tags= "Link,l=6")[1]);  
id_comb = combiner(id_mps, id_mpo; tags="id_comb");


I attached the output of versioninfo() as follows.

Julia Version 1.8.5
Commit 17cfb8e65e (2023-01-08 06:45 UTC)
Platform Info:
OS: Windows (x86_64-w64-mingw32)
CPU: 12 × AMD Ryzen 5 5600 6-Core Processor
LIBM: libopenlibm
LLVM: libLLVM-13.0.1 (ORCJIT, znver3)
Threads: 1 on 12 virtual cores

I appreciate your clear explanations which are easy to understand. Your answers to my beginner-level questions have been very helpful. Thank you for addressing my inquiries so well. ITensors seems to be a very impressive code.

There is an improved version of the 3S algorithm which is much faster and is actually easier to implement than the old 3S algorithm. For a new code, I would recommend implementing the new scheme (“post-expansion”) instead of 3S, as it has many advantages. It is easier to use because you don’t need to choose a mixing factor, instead just set some fraction of the basis that you dedicate to be expansion vectors (10% works well pretty much always). It is easier to implement because you can re-use the same contractions that you already use for the block updates and the Hamiltonian matrix-vector multiply. i.e. the contraction to find the vectors to add is just the usual Hamiltonian matrix-vector multiply, but with one of the environment tensors replaced by (weighted) random matrices. No paper yet, but there are some talk slides at that describe it, see the sections on ‘post-expansion’. The slides cover most of it except the weighting of the tensors, but this is probably not critical (I think in 3S, most people don’t use weighting anyway). If you need more details, email me and I can send you some notes, or wait for the paper which will be a couple of weeks (it would be out now, but I am traveling this month!).


Thanks a lot! A. I will implement it ^^. I was greatly inspired by your presentation at the recent conference SQAI. I heard that there’s an earthquake in Taiwan. Are you okay? Anyway, could i receive the notes? I could not find your email in google scholar ^^; My email is “” or “”. Thank you! Have a great trip!

@ianmcc is this the same procedure described in your arXiv posting [2403.00562] Comment on "Controlled Bond Expansion for Density Matrix Renormalization Group Ground State Search at Single-Site Costs" (Extended Version)?

It is not the main topic of the comment, but it is mentioned without any detail. The topic of the comment is an algorithm that is basically a faster version of 2-site DMRG. Post-expansion is basically the same ideas but applied to 3S instead: do the 3S mixing step in the null space, and use the RSVD instead of normal SVD.

Thanks for the clarification, sounds like these will be very useful techniques.