Hi,

I am trying to simulate a 1D-chain Spin 1/2 chain with a Hamiltonian containing 3-site terms, e.g. of the type S(i) S(i+1) S(i+2).

From my understanding the only way to implement this currently is to extend the basis per site to consist of 2 spins, spin1 and spin2, which allows me to write the 3-site term as a usual bond term.

I have implemented this in the attached file model_extended_basis.xml

Doing this then requires rewriting the 2-site terms, specifically the term S_plus(i) *S_minus(i+1)+h.c. will become a site term for even i S1_plus(i) *S2_minus(i) +h.c. and a bond term for uneven i S2_plus(i)*S1_minus(i+1) +h.c.

This implementation works flawlessly when using the full_diag routine. However, when using the dmrg routine, weird behaviour occurs.

Firstly, restricting the simulation to total S_z sectors works within the full_diag, but fails in dmrg for S_z !=0 .

Secondly, the bond term seems to be treated fine, i.e. only including the bond term in the Hamiltonian, the results agree with the full_diag results.

However, when including the site-term, the ground-state energy suddenly becomes 100 times the coupling constant, and does not agree with the expectation value of the bond-term anymore, i.e. setting all other terms in the Hamiltonian to 0 and measuring the expectation value of the bond-operator does not agree with E0.

I have attached the python scripts to start the simulation as dmrg_simulation_test.py and fulldiag_simulation_test.py

I would be very grateful if someone were able to help me.

Kind Regards