Dear all,
Now I’m studying some transport property based on multilayered system, with nonmagnetic layer/ferromagnetic layer (NM/FM) stack.
May I ask when creating a tight binding system, and setting hoppings by the following function:
syst[(lat(1, 0), lat(1, 1))] = t,
how to associate the hopping with a realistic system? Such as hoppings between FM and FM atoms, or between FM and NM atoms? I found in many literatures people just simply put t=1, but I assume the value should depend on the type of atoms.
My another question is, is it possible to define more than one type of hopping between two sites? For the tight-binding approximation of p electrons, there are two types of hopping integrals: Vppσ, and Vppπ, so I’m wondering if it’s possible to simulate the effect of both hoppings?
Thanks a lot for reading my questions! Cheers, Yu Li
Dear Yu,
In kwant you simulate models. So it is up to you to decide what corresponds best to your real system.
Yes you can change the value of the hopping depending on the justification you may find. It is also possible to put it as a parameter and do your study as a function of the value of "t". You can also choose it as a matrix in order to take into account the different orbitals in your system
I hope this helps, Adel
On Wed, Feb 5, 2020 at 11:22 AM Yu Li yu.li@postgrad.manchester.ac.uk wrote:
Dear all,
Now I’m studying some transport property based on multilayered system, with nonmagnetic layer/ferromagnetic layer (NM/FM) stack.
May I ask when creating a tight binding system, and setting hoppings by the following function:
syst[(lat(1, 0), lat(1, 1))] = t,
how to associate the hopping with a realistic system? Such as hoppings between FM and FM atoms, or between FM and NM atoms? I found in many literatures people just simply put t=1, but I assume the value should depend on the type of atoms.
My another question is, is it possible to define more than one type of hopping between two sites? For the tight-binding approximation of p electrons, there are two types of hopping integrals: Vppσ, and Vppπ, so I’m wondering if it’s possible to simulate the effect of both hoppings?
Thanks a lot for reading my questions! Cheers, Yu Li
Dear Yu Li, It seems that you deal within Slater Koster framework. I would like to share with you a enclosed example (I usual give to my student) which explains how to get the band structure of a real case . Let as assume that we have a flat 2d shape of silicene (buckling is not considered for simplicity) and we consider that we know the : Vppσ, as well as Vppπ. To sum up let us say (we have Es, px, py and pz like-orbitals) For the consideration of Fm and NM we have to generalize this example to your case. Hop this will help best, Adel
## ------ IMPORTING PACKAGES ------------------------- import numpy as np import tinyarray from numpy import pi,sqrt,arccos import kwant from kwant.wraparound import wraparound,plot_2d_bands from matplotlib import pyplot import scipy
### ----------------------------------------------------------------------------------- ### -------------------------- THE MAIN USEFUL FUNCTIONS------------------------------- ### ----------------------------------------------------------------------------------- def wraparound_syst(syst): return kwant.wraparound.wraparound(syst).finalized() ## or wraparound(syst).finalized()
def get_A(syst): B = np.asarray(syst._wrapped_symmetry.periods).T #print("============ np.linalg.pinv(B).T =============") #print(np.linalg.pinv(B).T) return np.linalg.pinv(B).T
def momentum_to_lattice(k, syst): """Transform momentum to the basis of reciprocal lattice vectors. See https://en.wikipedia.org/wiki/Reciprocal_lattice#Generalization_of_a_dual_la... """ A = get_A(syst) #print("============ A =============") #print(A) #print(scipy.linalg.lstsq(A, k)) k, residuals = scipy.linalg.lstsq(A, k)[:2] if np.any(abs(residuals) > 1e-7): raise RuntimeError("Requested momentum doesn't correspond" " to any lattice momentum.") #print("============ k =============") #print(k) return k
### for 2D, we omit k_z: k =(kx, ky). for 3D, we set k =(kx, ky, kz) def energies(k, syst): """ THE SOLUTION OF THE EIGENVALUE PROBLEM """ ## ----------- FOR 2D symmetry -------------- k_x, k_y = momentum_to_lattice(k, syst) params = {'k_x': k_x, 'k_y': k_y} ## ----------- FOR 3D symmetry -------------- ## we use #k_x, k_y, k_z = momentum_to_lattice(k, syst) #params = {'k_x': k_x, 'k_y': k_y, 'k_z': k_z} H = syst.hamiltonian_submatrix(params=params) energies = np.sort(np.linalg.eigvalsh(H)) return energies ## or np.array(energies)
### ----------------------------------------------------------------------------------- ### ------- THE TIGHT-BINDING PARAMETERS WITHIN SLATER FRAMEWORK ---------------------- ### ----------------------------------------------------------------------------------- ## --------------------- For example silicene (buckling is not considered for simplicity reason) Es_a = -4.0497; Es_b = -4.0497; ## --------------------- Ep_a = 1.0297; Ep_b = 1.0297; ## --------------------- Vss_sig_aa = -0.0000 Vss_sig_bb = -0.0000 Vss_sig_ab = -2.0662 ## --------------------- Vsp_sig_aa = +0.0000 Vsp_sig_bb = +0.0000 Vsp_sig_ab = +2.0850 ## -------------------- Vpp_sig_aa = +0.8900 Vpp_sig_bb = +0.8900 Vpp_sig_ab = +3.1837 ## -------------------- Vpp_pi_aa = -0.3612 Vpp_pi_bb = -0.3612 Vpp_pi_ab = -0.9488
## --------------------------------------------- ## ====== Onsites for the system =============== def Onsite_Mo(site): """ This function is used when we take into account masse term wher Ep_a = - Ep_b """ if site.family == a : ## ----------------------------------------------------------------- Epx, Epy, Epz = Ep_a*np.ones(3) Es = Es_a Onsite_Mat = tinyarray.array([ [ Es, 0, 0, 0], [ 0, Epx, 0, 0], ## Epx = Ep [ 0, 0, Epy, 0], ## Epy = Ep [ 0, 0, 0, Epz], ## Epz = Ep ]) else: ## ------------------------------------------------------------------ ## If we include the masse terme we set Epx, Epy, Epz = -1*Ep_b*np.ones(3) Epx, Epy, Epz = Ep_b*np.ones(3) Es = Es_b Onsite_Mat = tinyarray.array([ [ Es, 0, 0, 0], [ 0, Epx, 0, 0], ## Epx = Ep [ 0, 0, Epy, 0], ## Epy = Ep [ 0, 0, 0, Epz], ## Epz = Ep ])
return Onsite_Mat
## ------------------------------------------------ ## ====== Hopping for the system =============== def Hopp_Mo(site_1, site_2): ## ---------------------------------------------------- """ the hopping is between A-A B-B and A-B and given by slater koster parameters """ x_1, y_1 = site_1.pos[0], site_1.pos[1] x_2, y_2 = site_2.pos[0], site_2.pos[1] Norm_Vec = np.sqrt( (x_2 - x_1)**2 + (y_2 - y_1)**2) ## in cease of 3d we set #x_1, y_1, z_1 = site_1.pos[0], site_1.pos[1], site_1.pos[2] #x_2, y_2, z_2 = site_2.pos[0], site_2.pos[1], site_2.pos[2] #Norm_Vec = np.sqrt( (x_2 - x_1)**2 + (y_2 - y_1)**2+ (z_2 - z_1)**2)
## --------------------- ## Cosinus Directions l = (x_2 - x_1)/Norm_Vec m = (y_2 - y_1)/Norm_Vec ## For 3D we set : n = (z_2 - z_1)/Norm_Vec ## In the actual case wedeal with 2d so n=0 n = 0
## ----------------------------------------------------
if (site_1.family == a and site_2.family == b): ## --------------------------------------------------- #print("===============") #print(l) #print(m) #print(n) S_S = Vss_sig_ab S_px = l*Vsp_sig_ab S_py = m*Vsp_sig_ab S_pz = n*Vsp_sig_ab ## ---------------------------------------------------- px_px = (l**2)*Vpp_sig_ab + (1-l**2)*Vpp_pi_ab px_py = (l*m)*(Vpp_sig_ab - Vpp_pi_ab) px_pz = (l*n)*(Vpp_sig_ab - Vpp_pi_ab) ## ---------------------------------------------------- py_py = (m**2)*Vpp_sig_ab + (1-m**2)*Vpp_pi_ab py_pz = (m*n)*(Vpp_sig_ab - Vpp_pi_ab) ## ---------------------------------------------------- pz_pz = (n**2)*Vpp_sig_ab + (1-n**2)*Vpp_pi_ab ## ----------------------------------------------------
Hop_Mat = tinyarray.array([ [ S_S, S_px, S_py, S_pz], [-S_px, px_px, px_py, px_pz], ## -pya_Sb = -Sa_pxb [-S_py, px_py, py_py, py_pz], ## pya_pxb = pxa_pyb [-S_pz, px_pz, py_pz, pz_pz], ## pza_pxb = pxa_pzb ])
elif (site_1.family == a and site_2.family == a):
## --------------------------------------------------- S_S = Vss_sig_aa S_px = l*Vsp_sig_aa S_py = m*Vsp_sig_aa S_pz = n*Vsp_sig_aa ## ---------------------------------------------------- px_px = (l**2)*Vpp_sig_aa + (1-l**2)*Vpp_pi_aa px_py = (l*m)*(Vpp_sig_aa - Vpp_pi_aa) px_pz = (l*n)*(Vpp_sig_aa - Vpp_pi_aa) ## ---------------------------------------------------- py_py = (m**2)*Vpp_sig_aa + (1-m**2)*Vpp_pi_aa py_pz = (m*n)*(Vpp_sig_aa - Vpp_pi_aa) ## ---------------------------------------------------- pz_pz = (n**2)*Vpp_sig_aa + (1-n**2)*Vpp_pi_aa ## ----------------------------------------------------
Hop_Mat = tinyarray.array([ [ S_S, S_px, S_py, S_pz], [-S_px, px_px, px_py, px_pz], ## -pya_Sb = -Sa_pxb [-S_py, px_py, py_py, py_pz], ## pya_pxb = pxa_pyb [-S_pz, px_pz, py_pz, pz_pz], ## pza_pxb = pxa_pzb ]) elif (site_1.family == b and site_2.family == b):
## ---------------------------------------------------- S_S = Vss_sig_bb S_px = l*Vsp_sig_bb S_py = m*Vsp_sig_bb S_pz = n*Vsp_sig_bb ## ---------------------------------------------------- px_px = (l**2)*Vpp_sig_bb + (1-l**2)*Vpp_pi_bb px_py = (l*m)*(Vpp_sig_bb - Vpp_pi_bb) px_pz = (l*n)*(Vpp_sig_bb - Vpp_pi_bb) ## ---------------------------------------------------- py_py = (m**2)*Vpp_sig_bb + (1-m**2)*Vpp_pi_bb py_pz = (m*n)*(Vpp_sig_bb - Vpp_pi_bb) ## ---------------------------------------------------- pz_pz = (n**2)*Vpp_sig_bb + (1-n**2)*Vpp_pi_bb ## ----------------------------------------------------
Hop_Mat = tinyarray.array([ [ S_S, S_px, S_py, S_pz], [-S_px, px_px, px_py, px_pz], ## -pya_Sb = -Sa_pxb [-S_py, px_py, py_py, py_pz], ## pya_pxb = pxa_pyb [-S_pz, px_pz, py_pz, pz_pz], ## pza_pxb = pxa_pzb ])
return Hop_Mat
## -------------------------------- Make a lattice ------------------------------------- ## ===================================================================================== def make_2dHonycomb(): a_sisi = 1.42;
"""Create a builder for a periodic infinite sheet of graphene.""" a1 = a_sisi*sqrt(3) lat = kwant.lattice.general([[a1, 0], [a1/2, a1*sqrt(3)/2]], # lattice vectors [[ 0, 0], [ 0, a1/sqrt(3)]], # Positions name = ['a', 'b']) # sublatice'names a, b = lat.sublattices
## ---------- For unit cell ---------------------------------------------------------------- sym = kwant.TranslationalSymmetry(lat.vec((1, 0)), lat.vec((0, 1))) sys = kwant.Builder(sym) sys[a.shape(lambda p: True, (0, 0))] = Onsite_Mo sys[b.shape(lambda p: True, (0, 0))] = Onsite_Mo sys[a.neighbors(1)] = Hopp_Mo sys[b.neighbors(1)] = Hopp_Mo hoppings_ab = (((0, 0), a, b), ((0, +1), a, b), ((-1, +1), a, b)) sys[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings_ab]] = Hopp_Mo
#sys[lat.neighbors(2)] = ## -------------------------------------------------------------------------------------------
return sys, (a, b)
### ------------------- DEFINING THE UNIT CELL -------------- ### =========================================================
def hop_lw(site1, site2): return 0.01 if site1.family == site2.family else 0.03
def hopping_color(site1,site2): return 'r' if site1.family==site2.family else 'k' syst, (a, b) = make_2dHonycomb() kwant.plot(syst, site_size=0.1, site_lw=0.012, hop_color=hopping_color, hop_lw=hop_lw) pyplot.show()
### =============================== unit_cell = wraparound_syst(syst) ## =================================== High Symmetry Points========================================= ## ----- G ------- -------- K -------- ------ M ----- ------- G ------- ## x- 0 2*pi/(3*sqrt(3)*a_cc) 0 0 ## y- 0 2*pi/(3*a_cc) 2*pi/(3*a_cc) 0 ## ================================================================================================= a_sisi = 1.42; a_cc = a_sisi ac = sqrt(3)*a_cc d_GK = sqrt((2*pi/(3*sqrt(3)*a_cc)-(0))**2 + ( 2*pi/(3*a_cc)-0)**2) d_KM = sqrt(((0)-2*pi/(3*sqrt(3)*a_cc))**2 + ( 2*pi/(3*a_cc)-2*pi/(3*a_cc))**2) d_MG = sqrt(((0)-(0))**2 + ((0) - 2*pi/(3*a_cc))**2) print("the distance between G-K is " , d_GK ) print("the distance between K-M is " , d_KM ) print("the distance between M-G is ", d_MG ) ## ------------------------------------------------- ## we set the nb for np.linspace(, ,nb) for each path number = 400 ## we set this for the long path which is 1.7 for K1-M
nb_GK = number ### number*(d_GK/d_GK)= number## since number sinse GK is the longest path nb_KM = number*(d_KM/d_GK) nb_MG = number*(d_MG/d_GK) ## --------------------------------------------------------------------- kx_KG = np.linspace( 2*pi/(3*sqrt(3)*a_cc), 0, int(nb_GK) ) ky_KG = np.linspace( 2*pi/(3*a_cc), 0, int(nb_GK) ) ## --------------------------------------------------------------------- kx_GM = np.linspace( 0, 2*pi/(3*a_cc), int(nb_MG) ) ky_GM = np.linspace( 0, 2*pi/(3*a_cc), int(nb_MG) ) ## --------------------------------------------------------------------- kx_MK = np.linspace( 2*pi/(3*a_cc), 2*pi/(3*sqrt(3)*a_cc), int(nb_KM) ) ky_MK = 2*pi/(3*a_cc)*np.ones(int(nb_KM)) ## ---------------------------------------------------------------------
## For Kx Kx_path=[] Kx_path = np.append(Kx_path, kx_KG) Kx_path = np.append(Kx_path, kx_GM) Kx_path = np.append(Kx_path, kx_MK) ## For Kx Ky_path=[] Ky_path = np.append(Ky_path, ky_KG) Ky_path = np.append(Ky_path, ky_GM) Ky_path = np.append(Ky_path, ky_MK)
## -------------------------------------------------------------------------------- ## Extracting the energy bands along high symmetry point ========================== Energy_bands = [energies((kx, ky), unit_cell) for (kx, ky) in zip(Kx_path, Ky_path)]
## for HSP plot, we can just set a random vectors wich has the size of Energy_bands ## Here will will show G-K-M-K so we do not need the exact values of reciprocal vector size_HSP = np.size(np.array(Energy_bands)[:,0]) HSP_For_Plot = np.linspace(0, 1, size_HSP)
### --------------------------------------------------------------------------------------------- ###====================== SETTING THE PLOT PARAMETERS FOR G-K-M-G =============================== ### ----------------------- ---------------------------------------------------------------------- fig = pyplot.figure() #axes.axis('equal') axes = fig.add_subplot(1, 1, 1) for i in range(np.size(np.array(Energy_bands)[0])): axes.plot(HSP_For_Plot, np.array(Energy_bands)[:,i])
## ------------------------------------------------------- ## Setting the Horizontal lines for the high symmetry path axes.axvline(x = HSP_For_Plot[0], linewidth=2, color='k') axes.axvline(x = HSP_For_Plot[np.size(kx_KG)], linewidth=2, color='k') axes.axvline(x = HSP_For_Plot[np.size(kx_KG)+ np.size(kx_GM)], linewidth=2, color='k') axes.axvline(x = HSP_For_Plot[np.size(kx_KG)+ np.size(kx_GM)+ np.size(kx_MK)-1], linewidth=2, color='k')
axes.axhline(y = 0, color='k')
## ---------------------------------------------------------------------------- ## for HSP plot, we can just set a random vecto wich has the size of Energy_bands ## the normalized position for each high symmetry ponits are defined as: dK = HSP_For_Plot [0] dKG = HSP_For_Plot[np.size(kx_KG)] dKGM = HSP_For_Plot[np.size(kx_KG)+np.size(kx_GM)] dKGMK = HSP_For_Plot[np.size(kx_KG)+np.size(kx_GM)+np.size(kx_MK)-1] HighSymmetry_positions = [dK , dKG, dKGM, dKGMK] HighSymmetry_names= ('$K$',r'$\Gamma$', '$M$', '$K$') ## calling xticks and plotting horizontal lines related to G, K, M, and G from matplotlib.pyplot import xticks locs, labels = xticks() axes.xaxis.set_visible(True) axes.yaxis.set_visible(True) xticks(HighSymmetry_positions, HighSymmetry_names)
Le mer. 5 févr. 2020 à 19:56, Abbout Adel abbout.adel@gmail.com a écrit :
Dear Yu,
In kwant you simulate models. So it is up to you to decide what corresponds best to your real system.
Yes you can change the value of the hopping depending on the justification you may find. It is also possible to put it as a parameter and do your study as a function of the value of "t". You can also choose it as a matrix in order to take into account the different orbitals in your system
I hope this helps, Adel
On Wed, Feb 5, 2020 at 11:22 AM Yu Li yu.li@postgrad.manchester.ac.uk wrote:
Dear all,
Now I’m studying some transport property based on multilayered system, with nonmagnetic layer/ferromagnetic layer (NM/FM) stack.
May I ask when creating a tight binding system, and setting hoppings by the following function:
syst[(lat(1, 0), lat(1, 1))] = t,
how to associate the hopping with a realistic system? Such as hoppings between FM and FM atoms, or between FM and NM atoms? I found in many literatures people just simply put t=1, but I assume the value should depend on the type of atoms.
My another question is, is it possible to define more than one type of hopping between two sites? For the tight-binding approximation of p electrons, there are two types of hopping integrals: Vppσ, and Vppπ, so I’m wondering if it’s possible to simulate the effect of both hoppings?
Thanks a lot for reading my questions! Cheers, Yu Li
-- Abbout Adel
Le mer. 5 févr. 2020 à 09:22, Yu Li yu.li@postgrad.manchester.ac.uk a écrit :
Dear all,
Now I’m studying some transport property based on multilayered system, with nonmagnetic layer/ferromagnetic layer (NM/FM) stack.
May I ask when creating a tight binding system, and setting hoppings by the following function:
syst[(lat(1, 0), lat(1, 1))] = t,
how to associate the hopping with a realistic system? Such as hoppings between FM and FM atoms, or between FM and NM atoms? I found in many literatures people just simply put t=1, but I assume the value should depend on the type of atoms.
My another question is, is it possible to define more than one type of hopping between two sites? For the tight-binding approximation of p electrons, there are two types of hopping integrals: Vppσ, and Vppπ, so I’m wondering if it’s possible to simulate the effect of both hoppings?
Thanks a lot for reading my questions! Cheers, Yu Li