Dear Wilson, The greens function in kwant is calculated between the total sites connected to the lead only. In other words, it is a submatrix of the total green function you want to express. We sometimes call it: the surface greens function.
Please check my comments in the code below and the post [1] with the answer by Joseph in the kwant forum: I hope this helps, Adel [1] https://mail.python.org/archives/list/[email protected]/message/7BUBGEVE6Z4S4N7Y56M7DHDVYHLI3MH5/ ############################################################################## import kwant lat = kwant.lattice.square() def make_system(r=5, t=-1): def circle(pos): x, y = pos return x**2 + y**2 < r**2 syst = kwant.Builder() syst[lat.shape(circle, (0, 0))] = 0 syst[lat.neighbors()] = t syst.eradicate_dangling() def lead_shape(pos): return -4 < pos[0] < 4 lead = kwant.Builder( kwant.TranslationalSymmetry((0, 1))) lead[lat.shape(lead_shape, (0, 0))] = 0 lead[lat.neighbors()] = -t syst.attach_lead(lead) syst.attach_lead(lead.reversed()) return syst syst = make_system() fsyst = syst.finalized() kwant.plot(fsyst), #Get H: H = fsyst.hamiltonian_submatrix(sparse=False) #get the green's function matrix at energy E=-1.2 greens_function=kwant.greens_function(fsyst,energy=-1.2).data print('The shape of the greens function matrix is:' ,greens_function.shape) index0 = fsyst.lead_interfaces[0] # the index of sites to wich the lead number 0 is connected. index1 = fsyst.lead_interfaces[1] # the index of sites to wich the lead number 1 is connected. print('the total number of sites connected to a lead is: N=', len(index0)+len(index1)) print('system size:', H.shape) print('the greens function is NxN, which is diffferent from the size of the matrix H') def plot_color(site_index): if site_index in index0: return 'g' if site_index in index1: return 'y' else: return 'b' kwant.plot(fsyst,site_color=plot_color), print('If you want the total greens functions that includes the sites that are not connected to a lead') print('you need to construct the matrix manually by adding zero elements at the right place') print('Or you can do it by adding ficticious leads on the other sites. (check the kwant forum)') On Tue, May 28, 2024 at 3:23 PM <[email protected]> wrote: > Dear kwant community: > > I am investigating the `kwant.greens_function`。 I would like to get the > lead self-energies and add to the original system Hamiltonian and then > inverse it, thus obtaining the retarded green's function. I tested, but the > resulted green's function is not identical to the one kwant calculated. I > want to know the reason. Is it because what I describe in the title? > > Below is a minimal code example: > > ``` > lat = kwant.lattice.square() > > def make_system(r=20, t=-1): > > def circle(pos): > x, y = pos > return x**2 + y**2 < r**2 > > syst = kwant.Builder() > syst[lat.shape(circle, (0, 0))] = 0 > syst[lat.neighbors()] = t > syst.eradicate_dangling() > > def lead_shape(pos): > return -12 < pos[0] < 12 > > lead = kwant.Builder( kwant.TranslationalSymmetry((0, 1))) > > lead[lat.shape(lead_shape, (0, 0))] = 0 > lead[lat.neighbors()] = -t > syst.attach_lead(lead) > return syst > > syst = make_system() > fsyst = syst.finalized() > kwant.plot(fsyst) > > > greens_function=kwant.greens_function(fsyst) > Hc = fsyst.hamiltonian_submatrix(sparse=False) > > index = fsyst.lead_interfaces[0] > Hc[np.ix_(index, index)] += greens_function.lead_info[0] > > Gr = np.linalg.inv(-Hc) > > Gr[np.ix_(index, index)]-greens_function.submatrix(0, 0) # should be a > zero matrix > > ``` > > > Bests, > Wilson > -- Abbout Adel
