Dear Kwant community,
I would like to calculate the DOS and conductivity in 2D system with Zeeman and Rashba spin-orbit interaction in a square lattice. Regarding the chapter 2.3.1 in KWANT tutorial I defined the function make_system(). Then, according to the chapter 2.9 in the tutorial, I used the kwant.kpm.conductivity() method to calculate the longitudinal and transverse conductivity in the system. But I see that my results are wrong and I do not understand some things:
1. question: in chapter 2.3.1 we did not need to define the norbs parameter but I suppose that KWANT "knows" that we have here two orbitals per site (considering the spin) because we use the Pauli matrices. But, to obtain some results using function kwant.kpm.conductivity(), I see that the argument norbs=2 in the system deifinition is required (without it I get the error: 'NoneType' object cannot be interpreted as an integer). Thus, my quesion is: how does KWANT interpret my system if I add or not this argument (nobrs=2)? I think this could be the problem causing my wrong results.
2. question: To use the method kwant.kpm.conductivity() we don't need to attach the leads. Maybe it is obvious and trivial quesion, but I do not understand how we can calculate the conductivity (eg. induced by the external electric field) in a closed system? And what about the boundary conditions?
Below I present the whole code.
#--------------------------------------------------------------------------------------
#--------------------------------------------------------------------------------------
import kwant
from matplotlib import pyplot
import tinyarray
# we define the identity matrix and Pauli matrices
sigma_0 = tinyarray.array([[1,0], [0,1]])
sigma_x = tinyarray.array([[0,1], [1,0]])
sigma_y = tinyarray.array([[0,-1j], [1j,0]])
sigma_z = tinyarray.array([[1,0], [0,-1]])
def make_system(t=1.8, alpha=0.1, e_z=0.0018, W=30, L=30, a=1):
lat = kwant.lattice.square(norbs=2) # norbs - the number of orbitals per site
syst = kwant.Builder()
syst[(lat(x,y) for x in range(L) for y in range(W))] = 4*t*sigma_0 + e_z*sigma_z # Zeeman energy adds to the on-site term
syst[kwant.builder.HoppingKind((a,0), lat, lat)] = -t*sigma_0 - 1j*alpha*sigma_y/(2*a) # hoppings in x-direction; (a,0) means hopping form (i,j) to (i+1,j)
syst[kwant.builder.HoppingKind((0,a), lat, lat)] = -t*sigma_0 + 1j*alpha*sigma_x/(2*a) # hoppings in y-direction; (0,a) means hopping form (i,j) to (i,j+1)
return syst
def plot_dos_and_curves(dos, labels_to_data):
pyplot.figure()
pyplot.fill_between(dos[0], dos[1], label="DoS [a.u.]",
alpha=0.5, color='gray')
for label, (x, y) in labels_to_data:
pyplot.plot(x, y, label=label, linewidth=2)
pyplot.legend(loc=2, framealpha=0.5)
pyplot.xlabel("energy [t]")
pyplot.ylabel("$\sigma [e^2/h]$")
pyplot.show()
def main():
syst = make_system()
kwant.plot(syst); # check if the system looks as intended
syst = syst.finalized()
spectrum = kwant.kpm.SpectralDensity(syst, rng=0, num_vectors=30) # object that represents the density of states for the system
energies, densities = spectrum()
# component 'xx'
cond_xx = kwant.kpm.conductivity(syst, alpha='x', beta='x', rng=0)
# component 'xy'
cond_xy = kwant.kpm.conductivity(syst, alpha='x', beta='y', rng=0)
energies = cond_xx.energies
cond_array_xx = np.array([cond_xx(e, temperature=0.01) for e in energies])
cond_array_xy = np.array([cond_xy(e, temperature=0.01) for e in energies])
plot_dos_and_curves(
(spectrum.energies, spectrum.densities * 8), # why multiplied by 8?
[
(r'Longitudinal conductivity $\sigma_{xx} / 4$',
(energies, cond_array_xx.real / 4)),
(r'Hall conductivity $\sigma_{xy}$',
(energies, cond_array_xy.real))],
)
# standard Python construct to execute 'main' if the program is used as a script
if __name__ == '__main__':
main()
#--------------------------------------------------------------------------------------
#--------------------------------------------------------------------------------------
Regards,
Anna Krzyzewska