Error when using kroneckerproduct option for GroupedSite

[JackStevens]

Hi there,

I'm trying to implement my own model and run TEBD on it with non-nearest neighbour interactions. To do this, I've split my lattice into unit cells of 4 sites, which then are all nearest-neighbour (the sites are hardcore bosons, which I've created by simply removing the JW entries in the FermionSite class, I hope this is correct?)

As such, I need to apply a series of operators using the add_coupling method in the model initialisation. Here is the initialisation of the sites (and the lattice);

Python: Select all

site = [Hardcore_Boson(conserve='N'), Hardcore_Boson(conserve='N'), Hardcore_Boson(conserve='N'), Hardcore_Boson(conserve='N')]
site2 = GroupedSite(site)
lat = Lattice([U], site, bc='open', bc_MPS='finite')

where U is the number of unit cells in the lattice (each unit cell containing 4 hardcore boson sites).

I want to try and use add_coupling with operators defined by the kroneckerproduct of the 'C' and 'Cd' operators already defined in the site class. However, when I try to run

Python: Select all

C0 = site2.get_op('C0')
Cd1 = site2.get_op('Cd1')
C2 = site2.get_op('C2')
Cd3 = site2.get_op('Cd3')

site2.kroneckerproduct([C0, Cd1, C2, Cd3])

I get the error

Python: Select all

ValueError: incompatible LegCharge
self                                  | other  
LegPipe(shape (2, 2, 2, 2)->16,       |  +1    
    qconj (+1, +1, +1, +1)->+1;       | 0 [[0] 
    block numbers (2, 2, 2, 2)->5)    | 1  [1]]
 +1     |  +1     |  +1     |  +1     | 2      
0 [[0]  | 0 [[0]  | 0 [[0]  | 0 [[0]  |        
1  [1]] | 1  [1]] | 1  [1]] | 1  [1]] |        
2       | 2       | 2       | 2       |        
)  

I'm quite new to TeNPy, so I'm sure that I'm making a silly error, or asking it to do something that it isn't designed to do, but I can't seem to figure out what it is. Please excuse my being a total beginner, any help would be much appreciated.

Kind Regards,
Jack.

[Johannes]

Hi Jack, welcome to TeNPy!

(the sites are hardcore bosons, which I've created by simply removing the JW entries in the FermionSite class, I hope this is correct?)

The more straightforward way for HardcodeBosons would be to use the tenpy.networks.site.BosonSite with the default Nmax=1.

Here is the initialisation of the sites (and the lattice);

Are you using site or site2 for the lattice? Both of them are actually legit:

• Use site2 inside the lattice. This follows the implementation of the tenpy.models.spins_nnn.SpinChainNNN
• Use site inside the lattice, such that you have in total lat.N_sites == 4* U.
  Implement the couplings by setting the u=0, 1, 2, 3 explicitly when calling add_coupling_term. This would be more along the lines of tenpy.models.spins_nnn.SpinChainNNN2; in particular, if your model is actually translation invariant by one (ungrouped) site like the SpinChainNNN2 you don't need to distinguish the different u yet.
  So far, this generates a Model with couplings beyond nearest neighbors, but you can modify the __init__ to call group_sites, i.e. in your case self.group_sites(n=4). This should happen after the generation of the MPO, but before finishing the initialization of the tenpy.models.model.NearestNeighborModel.
  If you need help with that, let me know. (It could be helpful to add a model parameter for this option.)

That being said, lets come to your error:

However, when I try to run

Code: Select all

C0 = site2.get_op('C0')
Cd1 = site2.get_op('Cd1')
C2 = site2.get_op('C2')
Cd3 = site2.get_op('Cd3')

site2.kroneckerproduct([C0, Cd1, C2, Cd3])

The kroneckerproduct expects 4 operators, which are onsite operators on each of the 4 sites, not the grouped site. However, the C0, C1, C2, Cd3 operators are take from the grouped site and hence are already operators living on the grouped site, e.g. C0 represents \(c_0 \otimes1_1 \otimes 1_2 \otimes 1_3 \). Hence, the correct code would be:

Code: Select all

C0 = site[0].get_op('C')
Cd1 = site[1].get_op('Cd')
C2 = site[2].get_op('C2')
Cd3 = site[3].get_op('Cd3')
CCdCCd = site2.kroneckerproduct([C0, Cd1, C2, Cd3])

However, there is a much simpler alternative: simply use the actual product of the defined operators in site2:
\[( c_0 \otimes1_1 \otimes 1_2 \otimes 1_3 ) ( 1_0 \otimes c^\dagger_1\otimes 1_2 \otimes 1_3 ) (1_0 \otimes 1_1 \otimes c_2 \otimes 1_3 ) ( 1_0 \otimes1_1\otimes 1_2 \otimes c^\dagger_3 ) = ( c_0 \otimes c^\dagger_1\otimes c_2 \otimes c^\dagger_3 ) \]
Hence, you can get the operator you need with just a single line

Code: Select all

CCdCCd = site2.get_op('C0 Cd1 C2 Cd3')

This works because space-separated operator names are viewed as products of operators by get_op.

I hope this helps. Let us know if you need further assistance

[JackStevens]

Dear Johannes,

Thanks for your quick and detailed reply, I've implemented those changes and now have a much cleaner model defined. I'm going to go for the implementation of the SpinChainNNN example, and I'm sure I'll be back soon with questions about the implementation of the couplings. Thanks very much for your help.

Kind Regards,
Jack Stevens.

[JackStevens]

Dear Johannes,

As predicted, I'm having some trouble with the self.group_sites option, specifically where to call it in my model class. Below is the code that initialises my model

Python: Select all

class Iomin_Model(CouplingMPOModel, NearestNeighborModel):
    """Implementation of the 4-site Iomin Hamiltonian, suitable for TEBD
    The Hamiltonian is:
        
     \hat{H}=\sum_{k}\Tilde{h}\omega_k \hat{C}_{k}^\dagger\hat{C}_{k} +
     \frac{{\Tilde{h}}^2\beta}{2}\sum_{{\bf k}}A_{\bf k} \hat{C}_{k_1}^\dagger\hat{C}_{k_2}^\dagger\hat{C}_{k_3}\hat{C}_{k_4}
    
    Parameters
    -----------
    L : int
        Length of the chain.
        
    beta : float
        Non-linearity parameter
        
    omega_k : float (elements of array)
        Eigenenergies of the single-particle problem, from the mathematica notebook.
        
    A_k : float (elements of array)
        Amplitudes of eigenstates, again given by mathematica notebook.
        
        
    bc_MPS : {'finite'}
        MPS boundary conditions, finite only for our purposes.
    """
    def __init__(self, model_params):
        CouplingMPOModel.__init__(self, model_params)

        
    def init_sites(self, model_params):
        bosonsite = Hardcore_Boson_Site(conserve='N')
        site = [bosonsite, bosonsite, bosonsite, bosonsite]
        
        return site
       
        
    def init_lattice(self, model_params):
        sites = self.init_sites(model_params)
        #define the lattice
        lat = Lattice([U], sites, bc='open', bc_MPS='finite')
        
        return lat
    
    def init_terms(self, model_params):
        name = "4-site Iomin"
        L = get_parameter(model_params, 'L', 4, name)
        U = get_parameter(model_params, 'U', 1, name)
        beta = get_parameter(model_params, 'beta', 1., name)
        omega_k = get_parameter(model_params, 'omega_k', 1., name, asarray=True)
        A_k = get_parameter(model_params, 'A_k', 1., name, asarray=True)
        
        for i in range(0, L, 4):
            self.add_onsite(omega_k[i], i, 'N')
            self.add_onsite(omega_k[i+1], i+1, 'N')
            self.add_onsite(omega_k[i+2], i+2, 'N')
            self.add_onsite(omega_k[i+3], i+3, 'N')
            

        self.add_coupling(4*beta*A_k[1][3][0][2], 0, 'B Bd B Bd', 0, 'B Bd B Bd', 1)

   

I've added only one coupling, to test it. Without the self.group_sites option, I get the error

Python: Select all

  File "C:\Users\Jack-Surface\tenpy\networks\terms.py", line 632, in to_nn_bond_Arrays
    raise ValueError(msg.format(i=i, j=j2))

ValueError: Can't give nearest neighbor H_bond for long-range 0-4

Which is expected, as I've given it a coupling term that is not nearest neighbour. However, I can't seem to find the right place to call the group_sites attribute? Calling it before/after the CouplingMPOModel in __init__ doesn't seem to work. I've also tried calling it in the CouplingMPOModel class itself, after the initialisation of the MPO but before the NearestNeighborModel, but that also throws up some errors, namely

Python: Select all

  File "C:\Users\Jack-Surface\tenpy\models\model.py", line 202, in group_sites
    old_L = len(self.H_bond)

AttributeError: 'Iomin_Model' object has no attribute 'H_bond'

Any help would be much appreciated.

Kind Regards,
Jack Stevens.

[Johannes]

I've just added a the function from_MPOModel, which has an example how to use it.
Does that help?

[JackStevens]

Dear Johannes,

That is perfect, thank you! I now understand, and am able to implement my desired couplings. Thanks very much for your helpful advice.

Kind Regards,
Jack Stevens.