I am trying to find the ground state and calculate the ground state energy of a spinless fermi system using GroupedSite and DMRG. However, the problem is that I need the number of impurity to be fixed, while actually it always change and finally DMRG always give me the global minimum with another number of impurity and number of host particles. For example, when I set the number of host fermions and of the impurity to be 10 and 1 respectively in the initial state, DMRG would give me the ground state with 5.5 fermions and 5.5 impurities, no matter the impurity is fermion or boson. How could I keep the number of host particles and impurities conserved during DMRG?
Thanks.
The definition of the model is
Code: Select all
#define the model with grouped sites for fermionic impurity
class GroupedBosonImpModel(CouplingMPOModel):
def init_sites(self, model_params):
conserve_F = model_params.get('conserve_F', 'N')
filling_F = model_params.get('filling_F',0.5)
conserve_B = model_params.get('conserve_B', 'N')
filling_B = model_params.get('filling_B',0)
Nmax = model_params.get('Nmax',1)
model_label = model_params.get('label',None)
model_charges = model_params.get('charges','same')
HostSite = FermionSite(conserve_F,filling_F)
ImpSite = BosonSite(Nmax,conserve_B,filling_B)
site = GroupedSite([HostSite,ImpSite])
NH=site.N0
NI=site.N1
N_product=npc.tensordot(NH,NI,[1,0])
site.add_op('N_product',N_product)
return site
def init_terms(self, model_params):
#confirm that the input lattice is a chain
if not isinstance(self.lat,Chain):
raise ValueError("The input lattice is not a Chain")
#get the coupling constants, t*c_{i+1}^\dagger c_i,
# h*x_{i+1}^\dagger x_i
# V*n_e*n_e and g*n_x*n_e
# mue, mux the potential for particle number conservation
t = model_params.get('t', 1.)
h = model_params.get('h', 1.)
V = model_params.get('V', 0.)
g = model_params.get('g', 0.)
mue = model_params.get('mue', 1.)
mux = model_params.get('mux', 1.)
#implement chemical potential
self.add_onsite(-mue,0,'N0',plus_hc=True)
self.add_onsite(-mux,0,'N1',plus_hc=True)
#implement e-imp interaction
self.add_onsite(-g,0,'N_product',plus_hc=True)
#implement e-e interaction
self.add_coupling(-V,0,'dN0',0,'dN0',1,plus_hc=True)
#implement kinetic energy
self.add_coupling(-t,0,'Cd0',0,'C0',1,plus_hc=True)
self.add_coupling(-h,0,'Bd1',0,'B1',1,plus_hc=True)
The code for dmrg is
Code: Select all
p_state=[[2],[0]]
for i in range(1,N):
p_state=p_state+[[2],[0]]
import time
start_time = time.time()
psi_trial=MPS.from_lat_product_state(TestModel.lat,p_state)
Nimp_0=psi_trial.expectation_value('N1')
Ne_0=psi_trial.expectation_value('N0')
print('Nimp_0=',Nimp_0.sum())
print('Ne_0=',Ne_0.sum())
dmrg_params_groupedsite = {
'mixer': None, # setting this to True helps to escape local minima
'max_E_err': 1.e-10,
'trunc_params': {
'chi_max': 10,
'svd_min': 1.e-10,
},
'verbose': True,
'combine': True
}
dmrg_params = {
'mixer': None,
'max_E_err': 1.e-10,
}
eng_GS = dmrg.TwoSiteDMRGEngine(psi_trial, TestModel, dmrg_params_groupedsite)
E, psi = eng_GS.run() # the main work; modifies psi in place
print("--- %s seconds ---" % (time.time() - start_time))
print('E=', E)
Nimp=psi.expectation_value('N1')
Ne=psi.expectation_value('N0')
print('Nimp=',Nimp.sum())
print('Ne=',Ne.sum())
Nimp_0= 0.0
Ne_0= 10.0
final DMRG state not in canonical form up to norm_tol=1.00e-05: norm_err=2.08e-02
--- 10.078522682189941 seconds ---
E= -40.299434837756976
Nimp= 5.000000000000001
Ne= 5.000000000000001