Dear users and developers, I would like to attach Büttiker’s virtual lead (the leads with zero net currents) to all sites in the scattering region to mimic dephasing processes. The attached code works fine as long as the switch 'vl=False' (without attaching any virtual lead). I get equal conductance for up and down spin. As I understood, attaching virtual leads can be done via the following builder class in kwant: https://kwant-project.org/doc/1/reference/generated/kwant. builder.SelfEnergyLead#kwant.builder.SelfEnergyLead I have the retarded self energy of each virtual lead (- 1j*0.0025) that I can replace with '*selfenergy_func*' in this builder class (am I right?). Now if I change 'vl=True' : 1- I get an error message, related to symmetry but there is no argument to define symmetry... 2- how can I ensure zero net currents through these virtual leads? Thanks in advance for 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 gd = 0.005 #the gamma_d, dephasing parameter 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() 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)) #the section that cause the problem vl=False if vl: def mount_vlead(syst, vlead_interface, norbs): """Mounts virtual lead to interfaces provided. :sys: kwant.builder.Builder An unfinalized system to mount leads :vlead_interface: sequence of kwant.builder.Site Interface of lead :norb: integer Number of orbitals in system hamiltonian. """ dim = len(vlead_interface)*norbs sea = zeros((dim, dim), dtype=float) - 1j*0.5*gd def selfenergy_func(energy, args=()): return sea vlead = kwant.builder.SelfEnergyLead(selfenergy_func, vlead_interface) syst.leads.append(vlead) #attaching virtual lead to all sites vli = syst.sites() mount_vlead(syst, vli, norbs=no) 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 sym = kwant.TranslationalSymmetry([0, 0, -a]) dn_lead = kwant.Builder(sym, conservation_law=-sigma_z) 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 sym1 = kwant.TranslationalSymmetry([0, 0, a]) up_lead = kwant.Builder(sym1, conservation_law=-sigma_z) 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 = [] datau = [] datad = [] for ie in range(-320,520): energy = ie * 0.001 smatrix = kwant.smatrix(syst, energy=energy) energies.append(energy) Gu=smatrix.transmission((1, 0), 0) Gd=smatrix.transmission((1, 1), 0) datau.append(Gu) datad.append(Gd) fig = pyplot.figure() pyplot.plot(energies, datau, 'r--') pyplot.plot(energies, datad, 'b:') pyplot.legend(['Gu', 'Gd'], loc='upper left') pyplot.xlim([-0.4,0.65]) pyplot.ylim([-0.03,1.0]) pyplot.show()
Dear Patrik,
I have some remarks for you:
1) In the definition of your variable 'sea' you are adding a scalar to a
matrix. This will give you a result different from what you expect. Start
by correcting this: gd ----> gd *eye(dim)
2) Call your function mount_vlead after you define your system and leads
vli = syst.sites()
mount_vlead(syst, vli, norbs=no)
3) If you mount a virtual lead, you will get an error when calling the
scattering matrix:
AttributeError: 'SelfEnergyLead' object has no attribute 'modes'
because the leads you added have no modes.
you can ask for a submatrix of the scattering matrix by specify only
your real leads:
smatrix(sys,energy, in_lead=[0,1],out_lead=[0,1])
This also will not help and you will get the same error.
In you are case, you are adding the same lead (so the same self energy)
for all your sites. You can avoid this just by
not adding any virtual lead and instead add directly your self energy to
your potential Ex:
sys[site1]=pot+self_energy
By doing this, your Hamiltonian becomes not Hermitian but this does not
prevent you from calculating the scattering matrix:
smatrix(sys,..........., check_hermiticity=False)
Now, this does not help you in specifying the current flowing through your
imaginary leads and put it to zero.
The only way I see to achieve that, is to use 1D real leads and use a
potential on each one and find the combination which makes the current =0.
For this purpose, you need to use the conductance matrix
conductance_matrix() #check documentation.
I hope this helps.
Adel
On Thu, Jul 27, 2017 at 2:49 PM, Patrik Arvoy
Dear users and developers,
I would like to attach Büttiker’s virtual lead (the leads with zero net currents) to all sites in the scattering region to mimic dephasing processes.
The attached code works fine as long as the switch 'vl=False' (without attaching any virtual lead). I get equal conductance for up and down spin.
As I understood, attaching virtual leads can be done via the following builder class in kwant: https://kwant-project.org/doc/1/reference/generated/kwant.bu ilder.SelfEnergyLead#kwant.builder.SelfEnergyLead
I have the retarded self energy of each virtual lead (- 1j*0.0025) that I can replace with '*selfenergy_func*' in this builder class (am I right?). Now if I change 'vl=True' :
1- I get an error message, related to symmetry but there is no argument to define symmetry... 2- how can I ensure zero net currents through these virtual leads?
Thanks in advance for 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 gd = 0.005 #the gamma_d, dephasing parameter 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()
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))
#the section that cause the problem vl=False if vl: def mount_vlead(syst, vlead_interface, norbs): """Mounts virtual lead to interfaces provided.
:sys: kwant.builder.Builder An unfinalized system to mount leads :vlead_interface: sequence of kwant.builder.Site Interface of lead :norb: integer Number of orbitals in system hamiltonian. """ dim = len(vlead_interface)*norbs sea = zeros((dim, dim), dtype=float) - 1j*0.5*gd
def selfenergy_func(energy, args=()): return sea
vlead = kwant.builder.SelfEnergyLead(selfenergy_func, vlead_interface) syst.leads.append(vlead)
#attaching virtual lead to all sites vli = syst.sites() mount_vlead(syst, vli, norbs=no)
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
sym = kwant.TranslationalSymmetry([0, 0, -a]) dn_lead = kwant.Builder(sym, conservation_law=-sigma_z)
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
sym1 = kwant.TranslationalSymmetry([0, 0, a]) up_lead = kwant.Builder(sym1, conservation_law=-sigma_z) 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 = [] datau = [] datad = []
for ie in range(-320,520): energy = ie * 0.001 smatrix = kwant.smatrix(syst, energy=energy) energies.append(energy) Gu=smatrix.transmission((1, 0), 0) Gd=smatrix.transmission((1, 1), 0) datau.append(Gu) datad.append(Gd)
fig = pyplot.figure() pyplot.plot(energies, datau, 'r--') pyplot.plot(energies, datad, 'b:') pyplot.legend(['Gu', 'Gd'], loc='upper left') pyplot.xlim([-0.4,0.65]) pyplot.ylim([-0.03,1.0]) pyplot.show()
-- Abbout Adel
Dear Adel,
Thank you for the comments!
Let's assume for the simplicity that I want to add only one virtual lead
and use the conductance matrix approach.
So I have already two real leads then I add the third lead (going to be
virtual) similar to the other two, to some site.
I calculate the conductance_matrix() and the voltage on each lead by
setting the current on the first lead to 1, second lead to -1
and the virtual lead to 0.
cond = array([[smatrix.transmission(i, j) for j in range(3)] for i in
range(3)])
voltage = linalg.solve(cond, [1, -1, 0])
How do I calculate the conductance then?
Regards
Patrik
On 27 July 2017 at 17:13, Abbout Adel
Dear Patrik,
I have some remarks for you:
1) In the definition of your variable 'sea' you are adding a scalar to a matrix. This will give you a result different from what you expect. Start by correcting this: gd ----> gd *eye(dim)
2) Call your function mount_vlead after you define your system and leads vli = syst.sites() mount_vlead(syst, vli, norbs=no)
3) If you mount a virtual lead, you will get an error when calling the scattering matrix: AttributeError: 'SelfEnergyLead' object has no attribute 'modes' because the leads you added have no modes. you can ask for a submatrix of the scattering matrix by specify only your real leads: smatrix(sys,energy, in_lead=[0,1],out_lead=[0,1]) This also will not help and you will get the same error.
In you are case, you are adding the same lead (so the same self energy) for all your sites. You can avoid this just by not adding any virtual lead and instead add directly your self energy to your potential Ex: sys[site1]=pot+self_energy
By doing this, your Hamiltonian becomes not Hermitian but this does not prevent you from calculating the scattering matrix: smatrix(sys,..........., check_hermiticity=False)
Now, this does not help you in specifying the current flowing through your imaginary leads and put it to zero.
The only way I see to achieve that, is to use 1D real leads and use a potential on each one and find the combination which makes the current =0. For this purpose, you need to use the conductance matrix
conductance_matrix() #check documentation.
I hope this helps. Adel
On Thu, Jul 27, 2017 at 2:49 PM, Patrik Arvoy
wrote: Dear users and developers,
I would like to attach Büttiker’s virtual lead (the leads with zero net currents) to all sites in the scattering region to mimic dephasing processes.
The attached code works fine as long as the switch 'vl=False' (without attaching any virtual lead). I get equal conductance for up and down spin.
As I understood, attaching virtual leads can be done via the following builder class in kwant: https://kwant-project.org/doc/1/reference/generated/kwant.bu ilder.SelfEnergyLead#kwant.builder.SelfEnergyLead
I have the retarded self energy of each virtual lead (- 1j*0.0025) that I can replace with '*selfenergy_func*' in this builder class (am I right?). Now if I change 'vl=True' :
1- I get an error message, related to symmetry but there is no argument to define symmetry... 2- how can I ensure zero net currents through these virtual leads?
Thanks in advance for 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 gd = 0.005 #the gamma_d, dephasing parameter 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()
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))
#the section that cause the problem vl=False if vl: def mount_vlead(syst, vlead_interface, norbs): """Mounts virtual lead to interfaces provided.
:sys: kwant.builder.Builder An unfinalized system to mount leads :vlead_interface: sequence of kwant.builder.Site Interface of lead :norb: integer Number of orbitals in system hamiltonian. """ dim = len(vlead_interface)*norbs sea = zeros((dim, dim), dtype=float) - 1j*0.5*gd
def selfenergy_func(energy, args=()): return sea
vlead = kwant.builder.SelfEnergyLead(selfenergy_func, vlead_interface) syst.leads.append(vlead)
#attaching virtual lead to all sites vli = syst.sites() mount_vlead(syst, vli, norbs=no)
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
sym = kwant.TranslationalSymmetry([0, 0, -a]) dn_lead = kwant.Builder(sym, conservation_law=-sigma_z)
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
sym1 = kwant.TranslationalSymmetry([0, 0, a]) up_lead = kwant.Builder(sym1, conservation_law=-sigma_z) 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 = [] datau = [] datad = []
for ie in range(-320,520): energy = ie * 0.001 smatrix = kwant.smatrix(syst, energy=energy) energies.append(energy) Gu=smatrix.transmission((1, 0), 0) Gd=smatrix.transmission((1, 1), 0) datau.append(Gu) datad.append(Gd)
fig = pyplot.figure() pyplot.plot(energies, datau, 'r--') pyplot.plot(energies, datad, 'b:') pyplot.legend(['Gu', 'Gd'], loc='upper left') pyplot.xlim([-0.4,0.65]) pyplot.ylim([-0.03,1.0]) pyplot.show()
-- Abbout Adel
Dear Patrik,
The lead you added plays the role of a dephasing channel and there is no
current flowing through it.
If the current you are applying in lead0 and lead1 is small enough to
consider linear response theory, the conductance will read :
g= I/(v1-v0)
v1, v0 are elements of the vector voltage you calculated.
Ps: Be careful, the matrix G (or Cond in your case) is not invertible. You
need to fix v2=0 and then use the 2x2 submatrix of G to get your vector
voltage.
The generalization to N leads is straightforward.
I hope this helps.
Adel
On Mon, Jul 31, 2017 at 4:14 PM, Patrik Arvoy
Dear Adel,
Thank you for the comments!
Let's assume for the simplicity that I want to add only one virtual lead and use the conductance matrix approach. So I have already two real leads then I add the third lead (going to be virtual) similar to the other two, to some site.
I calculate the conductance_matrix() and the voltage on each lead by setting the current on the first lead to 1, second lead to -1 and the virtual lead to 0.
cond = array([[smatrix.transmission(i, j) for j in range(3)] for i in range(3)]) voltage = linalg.solve(cond, [1, -1, 0])
How do I calculate the conductance then?
Regards Patrik
On 27 July 2017 at 17:13, Abbout Adel
wrote: Dear Patrik,
I have some remarks for you:
1) In the definition of your variable 'sea' you are adding a scalar to a matrix. This will give you a result different from what you expect. Start by correcting this: gd ----> gd *eye(dim)
2) Call your function mount_vlead after you define your system and leads vli = syst.sites() mount_vlead(syst, vli, norbs=no)
3) If you mount a virtual lead, you will get an error when calling the scattering matrix: AttributeError: 'SelfEnergyLead' object has no attribute 'modes' because the leads you added have no modes. you can ask for a submatrix of the scattering matrix by specify only your real leads: smatrix(sys,energy, in_lead=[0,1],out_lead=[0,1]) This also will not help and you will get the same error.
In you are case, you are adding the same lead (so the same self energy) for all your sites. You can avoid this just by not adding any virtual lead and instead add directly your self energy to your potential Ex: sys[site1]=pot+self_energy
By doing this, your Hamiltonian becomes not Hermitian but this does not prevent you from calculating the scattering matrix: smatrix(sys,..........., check_hermiticity=False)
Now, this does not help you in specifying the current flowing through your imaginary leads and put it to zero.
The only way I see to achieve that, is to use 1D real leads and use a potential on each one and find the combination which makes the current =0. For this purpose, you need to use the conductance matrix
conductance_matrix() #check documentation.
I hope this helps. Adel
On Thu, Jul 27, 2017 at 2:49 PM, Patrik Arvoy
wrote: Dear users and developers,
I would like to attach Büttiker’s virtual lead (the leads with zero net currents) to all sites in the scattering region to mimic dephasing processes.
The attached code works fine as long as the switch 'vl=False' (without attaching any virtual lead). I get equal conductance for up and down spin.
As I understood, attaching virtual leads can be done via the following builder class in kwant: https://kwant-project.org/doc/1/reference/generated/kwant.bu ilder.SelfEnergyLead#kwant.builder.SelfEnergyLead
I have the retarded self energy of each virtual lead (- 1j*0.0025) that I can replace with '*selfenergy_func*' in this builder class (am I right?). Now if I change 'vl=True' :
1- I get an error message, related to symmetry but there is no argument to define symmetry... 2- how can I ensure zero net currents through these virtual leads?
Thanks in advance for 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 gd = 0.005 #the gamma_d, dephasing parameter 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()
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))
#the section that cause the problem vl=False if vl: def mount_vlead(syst, vlead_interface, norbs): """Mounts virtual lead to interfaces provided.
:sys: kwant.builder.Builder An unfinalized system to mount leads :vlead_interface: sequence of kwant.builder.Site Interface of lead :norb: integer Number of orbitals in system hamiltonian. """ dim = len(vlead_interface)*norbs sea = zeros((dim, dim), dtype=float) - 1j*0.5*gd
def selfenergy_func(energy, args=()): return sea
vlead = kwant.builder.SelfEnergyLead(selfenergy_func, vlead_interface) syst.leads.append(vlead)
#attaching virtual lead to all sites vli = syst.sites() mount_vlead(syst, vli, norbs=no)
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
sym = kwant.TranslationalSymmetry([0, 0, -a]) dn_lead = kwant.Builder(sym, conservation_law=-sigma_z)
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
sym1 = kwant.TranslationalSymmetry([0, 0, a]) up_lead = kwant.Builder(sym1, conservation_law=-sigma_z) 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 = [] datau = [] datad = []
for ie in range(-320,520): energy = ie * 0.001 smatrix = kwant.smatrix(syst, energy=energy) energies.append(energy) Gu=smatrix.transmission((1, 0), 0) Gd=smatrix.transmission((1, 1), 0) datau.append(Gu) datad.append(Gd)
fig = pyplot.figure() pyplot.plot(energies, datau, 'r--') pyplot.plot(energies, datad, 'b:') pyplot.legend(['Gu', 'Gd'], loc='upper left') pyplot.xlim([-0.4,0.65]) pyplot.ylim([-0.03,1.0]) pyplot.show()
-- Abbout Adel
-- Abbout Adel
participants (2)
-
Abbout Adel
-
Patrik Arvoy