I started using TeNPy just a few days ago, and I appreciate all the work being done on it. Please forgive me for the following long, beginner questions, but I wasn't sure where else to go to for help.
1) I want to use DMRG to find the magnetic ground state on a triangular lattice for a spin-1/2 model with nearest-neighbor and next-nearest neighbor interactions (as seen in this Wikipedia page). I wrote the following code, but am not sure whether it is correct (when visualized, the lattice looks reasonable, but I am not sure about what other sanity checks I can perform). Do you have any advice on what I could do? Also, notes on general TeNPy-specific good practices are welcome (since I am still getting used to TeNPy).
Code: Select all
import numpy as np
import matplotlib.pyplot as plt
from tenpy.models.spins import SpinModel
from tenpy.networks.mps import MPS
from tenpy.algorithms import dmrg
from tenpy.models.lattice import Triangular
from tenpy.networks.site import SpinHalfSite
from tenpy.models.model import CouplingMPOModel
bc_xx = 'periodic'; bc_yy = 'periodic'; bc_MPSS = 'infinite'; # boundary conditions
coupling_axis = 'Sz' # adds exchange interaction along Sz axis
lx = 3; ly=3; tot_sites = lx*ly; # lattice
initial_state = ['up'] * (tot_sites) # initial state
j1= 1; j2 = 0.5; # exchange interactions
class MyModel(CouplingMPOModel):
def init_sites(self, model_params):
conserve = model_params.get('conserve', None)
site = SpinHalfSite(conserve=conserve)
return site
def init_lattice(self, model_params):
Lx = model_params.get('Lx', 3)
Ly = model_params.get('Ly', 3)
bc_x = model_params.get('bc_x', bc_xx)
bc_y = model_params.get('bc_y', bc_yy)
bc_MPS = model_params.get('bc_MPS', bc_MPSS)
lattice = Triangular(Lx, Ly, site=self.init_sites(model_params), bc=[bc_x, bc_y], bc_MPS=bc_MPS)
return lattice
def init_terms(self, model_params):
J1 = model_params.get('J1', 1.0)
J2 = model_params.get('J2', 0.0)
for u1, u2, dx in self.lat.pairs['nearest_neighbors']:
self.add_coupling(J1, u1, coupling_axis, u2, coupling_axis, dx)
for u1, u2, dx in self.lat.pairs['next_nearest_neighbors']:
self.add_coupling(J2, u1, coupling_axis, u2, coupling_axis, dx)
dmrg_params = {
'mixer': True,
'trunc_params': {
'chi_max': 10000,
'svd_min': 1.e-10,
},
'max_E_err': 1.e-10,
'verbose': 1,
}
model_params = {
'conserve': None,
'Lx': lx,
'Ly': ly,
'J1': j1,
'J2': j2,
'bc_MPS': bc_MPSS,
'bc_x': bc_xx,
'bc_y': bc_yy
}
M = MyModel(model_params)
lattice = M.init_lattice(model_params)
psi = MPS.from_product_state(lattice.mps_sites(), initial_state, bc=lattice.bc_MPS)
eng = dmrg.TwoSiteDMRGEngine(psi, M, dmrg_params)
info = eng.run() # the main work; modifies psi in place
E = info[0]; psi = info[1];
X = psi.expectation_value("Sx"); Y = psi.expectation_value("Sy"); Z = psi.expectation_value("Sz");
2) I want to test convergence of the ground state with respect to the number of states. However, almost all simulations I run (even up to 30 x 30 size, for both finite, semi-finite and infinite systems), use at most a maximum of chi = 6. Am I right in saying that this is the number of states used (bond dimension)? If so, is there a way for me to fix the number of states used, for me to do this type of convergence test? That is, I want to calculate the ground state energy for each number of states. Although I set chi_max to 10000 above, it doesn't go as anywhere close. I also tried chi_list, but to no avail. Any advice? Or just generally sanity checks I could use to make sure I am starting with a good ground state?
3) Finally, I want to play around and see whether I can use these simulations to check the type of magnetic ground state I have. For example, stripe, 120 AFM, FM, etc. Such as those in the image below (source):
![Image](https://www.researchgate.net/publication/342732797/figure/fig3/AS:921969021161472@1596825952331/Color-online-Top-view-of-the-effective-triangular-lattice-formed-by-V-ions-in-VXY.png)
Can anyone provide me a direction on how I could go about using TeNPy to classify simulation systems into those states? I naively thought that I could just calculate the expectation values of the spin operators at each site as I do with X, Y and Z in my code above. And then plot those as vectors (X,Y,Z) to visualize the spin. But, for the range of values I tested, I got only non-zero Z (and only values of 0.5 or -0.5, when I expected at least some cases to be 120 AFM, etc. Am I right in saying that it's more complicated than just getting the spin operator expectation values? If so, I would appreciate any directions I could use to start exploring. This paper uses TeNPy to achieve a similar end result. However, they use a so-called spin structure factor and momentum-space considerations to do so. I don't quite understand it yet (I am thinking I should perform a Fourier transform). But I would like to know what exactly I must do to intuitively classify these ground states.
I apologize for the long post, but I sincerely appreciate your time. Thanks again for all the work that goes into development!