Dear kwant users and developers, I am trying to calculate the up and down conductance in the following system, separately, to calculate spin polarization. I followed some of the steps suggested in the mailing list before (wavefuncation or current operator). I know for this system I should get the same conductance for up and down spins. I first calculate the total conductance But then 1- The wave function approach does not give the the up/down conductance similar to total. 2- the current-operator approach gives the error 'kwant.graph.core.EdgeDoesNotExistError'. I have attached the code in the following. If one runs the code it plot the system with leads correctly the first conductance plot is what one expects (if trans:) then the second plot is incorrect (if strans1:) and then gives the error (if strans:). I appreciate any help Patrik ----------------------------------------------------------------- from mpl_toolkits.mplot3d import Axes3D from scipy.spatial import * from matplotlib import rcParams from numpy import * from numpy.linalg import * import pickle import sys import os import string import heapq import kwant import tinyarray from matplotlib import pyplot chiral=True if chiral: p = pi/5 #phi t = 0.66 #theta a = 0.34 x = 1.4 e1 = 0 e2 = 0.3 t2=0.1 t1=-x*t2 t0 = 2 lam=-0.08 t_so1 = 0.01 #spin-orbit coupling param t_so2 = x*t_so1 #spin-orbit coupling param tl=tr=0.3 N = 30 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]]) no=2 #number of orbitals def sigma_v1(ap): # pauli metrix along the vertical axis value=sigma_z*cos(t)+sin(t)*(sigma_x*sin(ap)-sigma_y*cos(ap)) return value def sigma_v2(ap): # pauli metrix along the vertical axis value=sigma_z*cos(t)-sin(t)*(sigma_x*sin(ap)-sigma_y*cos(ap)) return value def family_color(sites): return 'black' #if site.family == sites def hopping_lw(site1, site2): return 0.08 class Amorphous(kwant.builder.SiteFamily): def __init__(self, coords): self.coords = coords super(Amorphous, self).__init__("amorphous", "",no) def normalize_tag(self, tag): try: tag = int(tag[0]) except: raise KeyError if 0 <= tag < len(coords): return tag else: raise KeyError def pos(self, tag): return self.coords[tag] coords=[(0.0000000000, 0.0000000000, 0.0000000000), (-0.1336881039, 0.4114496766, 0.3400000000), (-0.4836881039, 0.6657395614, 0.6800000000), (-0.9163118961, 0.6657395614, 1.0200000000), (-1.2663118961, 0.4114496766, 1.3600000000), (-1.4000000000, 0.0000000000, 1.7000000000), (-1.2663118961, -0.4114496766, 2.0400000000), (-0.9163118961, -0.6657395614, 2.3800000000), (-0.4836881039, -0.6657395614, 2.7200000000), (-0.1336881039, -0.4114496766, 3.0600000000), (0.0000000000, -0.0000000000, 3.4000000000), (-0.1336881039, 0.4114496766, 3.7400000000), (-0.4836881039, 0.6657395614, 4.0800000000), (-0.9163118961, 0.6657395614, 4.4200000000), (-1.2663118961, 0.4114496766, 4.7600000000), (-1.4000000000, 0.0000000000, 5.1000000000), (-1.2663118961, -0.4114496766, 5.4400000000), (-0.9163118961, -0.6657395614, 5.7800000000), (-0.4836881039, -0.6657395614, 6.1200000000), (-0.1336881039, -0.4114496766, 6.4600000000), (0.0000000000, -0.0000000000, 6.8000000000), (-0.1336881039, 0.4114496766, 7.1400000000), (-0.4836881039, 0.6657395614, 7.4800000000), (-0.9163118961, 0.6657395614, 7.8200000000), (-1.2663118961, 0.4114496766, 8.1600000000), (-1.4000000000, 0.0000000000, 8.5000000000), (-1.2663118961, -0.4114496766, 8.8400000000), (-0.9163118961, -0.6657395614, 9.1800000000), (-0.4836881039, -0.6657395614, 9.5200000000), (-0.1336881039, -0.4114496766, 9.8600000000), (-1.4000000000, 0.0000000000, 0.0000000000), (-1.2663118961, -0.4114496766, 0.3400000000), (-0.9163118961, -0.6657395614, 0.6800000000), (-0.4836881039, -0.6657395614, 1.0200000000), (-0.1336881039, -0.4114496766, 1.3600000000), (0.0000000000, -0.0000000000, 1.7000000000), (-0.1336881039, 0.4114496766, 2.0400000000), (-0.4836881039, 0.6657395614, 2.3800000000), (-0.9163118961, 0.6657395614, 2.7200000000), (-1.2663118961, 0.4114496766, 3.0600000000), (-1.4000000000, 0.0000000000, 3.4000000000), (-1.2663118961, -0.4114496766, 3.7400000000), (-0.9163118961, -0.6657395614, 4.0800000000), (-0.4836881039, -0.6657395614, 4.4200000000), (-0.1336881039, -0.4114496766, 4.7600000000), (0.0000000000, -0.0000000000, 5.1000000000), (-0.1336881039, 0.4114496766, 5.4400000000), (-0.4836881039, 0.6657395614, 5.7800000000), (-0.9163118961, 0.6657395614, 6.1200000000), (-1.2663118961, 0.4114496766, 6.4600000000), (-1.4000000000, 0.0000000000, 6.8000000000), (-1.2663118961, -0.4114496766, 7.1400000000), (-0.9163118961, -0.6657395614, 7.4800000000), (-0.4836881039, -0.6657395614, 7.8200000000), (-0.1336881039, -0.4114496766, 8.1600000000), (0.0000000000, -0.0000000000, 8.5000000000), (-0.1336881039, 0.4114496766, 8.8400000000), (-0.4836881039, 0.6657395614, 9.1800000000), (-0.9163118961, 0.6657395614, 9.5200000000), (-1.2663118961, 0.4114496766, 9.8600000000)] amorphous_lat = Amorphous(coords) syst = kwant.Builder() #adding the onsite and hopping to the system for i in range(N): syst[amorphous_lat(i)] = e1*sigma_0 syst[amorphous_lat(N+i)] = e2*sigma_0 syst[amorphous_lat(i), amorphous_lat(N+i)] = lam*sigma_0 if i > 0: syst[amorphous_lat(i), amorphous_lat(i-1)] = t1*sigma_0 + 1j*t_so1*(sigma_v1(i*p)+sigma_v1((i-1)*p)) syst[amorphous_lat(N+i),amorphous_lat(N+i-1)] = t2*sigma_0 + 1j*t_so2*(sigma_v2(i*p)+sigma_v2((i-1)*p)) # If we want to attach to vertical 1D chains to the system # we first add a site of the down lead to the scattering region lat=kwant.lattice.cubic(a, norbs=no) syst[lat(0, 0, -1)] = e1*sigma_0 syst[amorphous_lat(0), lat(0, 0, -1)] = tl*sigma_0 # We make a regular down lead and attach it to the system dn_lead = kwant.Builder(kwant.TranslationalSymmetry((0, 0, -a))) dn_lead[lat(0, 0, -2)] = e1*sigma_0 dn_lead[lat.neighbors()] = t0*sigma_0 syst.attach_lead(dn_lead) prim_vecs=tinyarray.array([(a,0.,0.),(0.,a,0.),(0.,0.,a)]) offset=tinyarray.array((-1.2663118961, 0.4114496766,0.0)) lat2=kwant.lattice.Monatomic(prim_vecs, offset, norbs=no) syst[lat2(0, 0, N)] = e1*sigma_0 syst[amorphous_lat(2*N-1), lat2(0, 0, N)] = tr*sigma_0 up_lead = kwant.Builder(kwant.TranslationalSymmetry((0, 0, a))) up_lead[lat2(0, 0, N+1)] = e1*sigma_0 up_lead[lat2.neighbors()] = t0*sigma_0 syst.attach_lead(up_lead) system=kwant.plot(syst, site_lw=0.1, site_color=family_color, hop_lw=hopping_lw) trans=True if trans: syst = syst.finalized() energies = [] data = [] for ie in range(-320,520): energy = ie * 0.001 smatrix = kwant.smatrix(syst, energy=energy) energies.append(energy) data.append(0.5*smatrix.transmission(1, 0)) fig = pyplot.figure(figsize=(6,2)) pyplot.plot(energies, data) pyplot.xlim([-0.32,0.52]) pyplot.ylim([-0.03,1.25]) pyplot.xlabel("energy [eV]") pyplot.ylabel("conductance [e^2/h]") pyplot.show() strans1=True if strans1: #syst = syst.finalized() energies = [] data = [] def oscle(ene, lead_nr): wfs=kwant.wave_function(syst, ene, check_hermiticity=True)(lead_nr) spin_current_z = 0 for psi in wfs: psi_start = psi[0 : 2] psi_end = psi[2 * 61: 2 * 61 + 2] #spin_current_z += -2 * imag(psi_end.conjugate().dot(sigma_z).dot(psi_start)) spin_current_z += abs(imag(psi_end.conjugate().dot(sigma_z).dot(psi_start))) return spin_current_z for ie in range(-320,520): energy = ie * 0.001 energies.append(energy) data.append(oscle(ene=energy, lead_nr=1)) fig = pyplot.figure(figsize=(6,2)) pyplot.plot(energies, data) pyplot.show() strans=True if strans: #syst = syst.finalized() J_spin = kwant.operator.Current(syst, sigma_z, where=[(amorphous_lat(0), amorphous_lat(N-2))], sum=True) all_wfs = kwant.wave_function(syst, energy=0.25)(1) spin_current_list = sum(J_spin(wf) for wf in all_wfs) print(spin_current_list)
Hi Patrik,
1- The wave function approach does not give the the up/down conductance similar to total.
I see that you do the following:
psi_start = psi[0 : 2] psi_end = psi[2 * 61: 2 * 61 + 2]
Do you know if the first and last sites sites as they are ordered in the finalized system are indeed the sites that you want? In your script you never refer to the 'sites' attribute of the finalized system, so I would guess that these are not the sites that you want; sites number 0 and 61 are probably not even connected.
2- the current-operator approach gives the error 'kwant.graph.core.EdgeDoesNotExistError'.
This means that you are trying to calculate the current between 2 sites that are not connected by a hopping.
J_spin = kwant.operator.Current(syst, sigma_z, where=[(amorphous_lat(0), amorphous_lat(N-2))], sum=True)
I see that you are trying to calculate the current between 'amorphous_lat(0)' and 'amorphous_lat(N-2)', and if I colour these sites when plotting the system I see indeed that they are at opposite ends of your double helix structure, and that they are not connected by a hopping. Happy Kwanting, Joe
Hi Joe,
Thank you for the reply!
As you mentioned they are not connected via a single hopping. I was trying
to get up/down conductance through the system
from one end to the other, similar to the conductance from lead 0 to lead 1.
What would be the equivalent of conductance from lead 0 to lead 1 (in the
'if trans'), but for either up or down spin?
Regards
Patrik
On 21 July 2017 at 12:12, Joseph Weston
Hi Patrik,
1- The wave function approach does not give the the up/down conductance similar to total.
I see that you do the following:
psi_start = psi[0 : 2] psi_end = psi[2 * 61: 2 * 61 + 2]
Do you know if the first and last sites sites as they are ordered in the finalized system are indeed the sites that you want? In your script you never refer to the 'sites' attribute of the finalized system, so I would guess that these are not the sites that you want; sites number 0 and 61 are probably not even connected.
2- the current-operator approach gives the error 'kwant.graph.core.EdgeDoesNotExistError'.
This means that you are trying to calculate the current between 2 sites that are not connected by a hopping.
J_spin = kwant.operator.Current(syst, sigma_z, where=[(amorphous_lat(0), amorphous_lat(N-2))], sum=True)
I see that you are trying to calculate the current between 'amorphous_lat(0)' and 'amorphous_lat(N-2)', and if I colour these sites when plotting the system I see indeed that they are at opposite ends of your double helix structure, and that they are not connected by a hopping.
Happy Kwanting,
Joe
Hi Patrik,
I was trying to get up/down conductance through the system from one end to the other, similar to the conductance from lead 0 to lead 1.
What would be the equivalent of conductance from lead 0 to lead 1 (in the 'if trans'), but for either up or down spin?
Ah, this is a new feature introduced in Kwant 1.3! The way to do this is to first tell Kwant that your leads have a conservation law associated with them, and to tell Kwant the basis that you want to use for the modes in the lead. Without this information the modes Kwant calculates in the leads will be an arbitrary superposition of spin up and spin down, and you will not easily be able to separate these. You provide this information with the 'conservation_law' parameter when creating the Builder for your leads. This parameter should be a Hermitian matrix with integer eigenvalues, who's eigenvectors you want to take as the basis for the modes. The modes will be ordered by their eigenvalues. This is explained in the docs[1]; the example given is for superconductivity, but the principle is the same. For your case a reasonable choice would be sigma_z: sym = kwant.TranslationalSymmetry([0, 0, -a]) dn_lead = kwant.Builder(sym, conservation_law=-sigma_z) … There is a minus sign in front of the sigma_z so that the "up" modes come first, and then the "down" modes. Then, when computing the transmission, instead of saying smatrix.transmission(1, 0) which will give you the total transmission from lead 0 to lead 1, you can say smatrix.transmission((1, 0), (0, 0)) to get the transmission from the spin up block (block 0) of lead 0 to the spin up block of lead 1. The first number in each pair is the lead number, and the second is the block index with respect to the conservation law you defined when creating the Builder (0 being spin up, and 1 being spin down in this case). Similarly you can use smatrix.transmission((1, 1), (0, 1)) for transmission from spin down to spin down, or any other combination to calculate the transmissions between different spin blocks. Hope that helps, Joe [1]: https://kwant-project.org/doc/1/tutorial/superconductors
Hi Joe,
Perfect, it works now.
Thanks a lot!
Regards
Patrik
On 21 July 2017 at 15:24, Joseph Weston
Hi Patrik,
I was trying to get up/down conductance through the system from one end to the other, similar to the conductance from lead 0 to lead 1.
What would be the equivalent of conductance from lead 0 to lead 1 (in the 'if trans'), but for either up or down spin?
Ah, this is a new feature introduced in Kwant 1.3! The way to do this is to first tell Kwant that your leads have a conservation law associated with them, and to tell Kwant the basis that you want to use for the modes in the lead. Without this information the modes Kwant calculates in the leads will be an arbitrary superposition of spin up and spin down, and you will not easily be able to separate these.
You provide this information with the 'conservation_law' parameter when creating the Builder for your leads. This parameter should be a Hermitian matrix with integer eigenvalues, who's eigenvectors you want to take as the basis for the modes. The modes will be ordered by their eigenvalues. This is explained in the docs[1]; the example given is for superconductivity, but the principle is the same. For your case a reasonable choice would be sigma_z:
sym = kwant.TranslationalSymmetry([0, 0, -a]) dn_lead = kwant.Builder(sym, conservation_law=-sigma_z) …
There is a minus sign in front of the sigma_z so that the "up" modes come first, and then the "down" modes.
Then, when computing the transmission, instead of saying
smatrix.transmission(1, 0)
which will give you the total transmission from lead 0 to lead 1, you can say
smatrix.transmission((1, 0), (0, 0))
to get the transmission from the spin up block (block 0) of lead 0 to the spin up block of lead 1. The first number in each pair is the lead number, and the second is the block index with respect to the conservation law you defined when creating the Builder (0 being spin up, and 1 being spin down in this case).
Similarly you can use
smatrix.transmission((1, 1), (0, 1))
for transmission from spin down to spin down, or any other combination to calculate the transmissions between different spin blocks.
Hope that helps,
Joe
[1]: https://kwant-project.org/doc/1/tutorial/superconductors
participants (2)
-
Joseph Weston
-
Patrik Arvoy