Dear Anton,
I have adopted my conservation law to my system (graphene/TMDs like
heterostructure). Band calculation is with correct agreement compared to
dft calculation.
Everything is fine unless we add hopping terms between spin up and spin
down. In this case, we cannot call spin transmission for spin_dn_dn; up_dn
and dn_up.
To answer my doubts, I have used graphene with a rashba spin-orbit coupling
(we are hopping elements between spin up and down like in my case) and the
same thing happens to me when calling
*dn_dn_trans = smatrix.transmission((1, 1), (0, 1)). *
It works only if we set params.rashba to zero.
*It is not correct to use the conservation law when we have coupling
between spins, like in rashba soc!*
I have attached below the script code of graphene with intrinsic and rashba
spin orbit terms
Thanks in advance
import kwant
import numpy as np
from math import sqrt
import matplotlib.pyplot as plt
from types import SimpleNamespace
import warnings
warnings.filterwarnings('ignore')
##
==========================================================================================================================
### ============ Pauli matrices ===============
### --------------------------------------------
sigma_0 = np.array([[1, 0], [0, 1]]) # identity
sigma_x = np.array([[0, 1], [1, 0]])
sigma_y = np.array([[0, -1j], [1j, 0]])
sigma_z = np.array([[1, 0], [0, -1]])
### --------- Onsite-------
def onsite(site, params):
return (params.Ep+params.mass)*sigma_0 if site.family ==a else
(params.Ep-params.mass)*sigma_0
### ============ Hopping terms for the diffent parts of the system
===============
###
---------------------------------------------------------------------------
def kinetic(site_1, site_2, params):
""" use it for nearest interaction"""
return params.hop*sigma_0
###
-------------------------------------------------------------------------
def Rashba_SOC(site_1, site_2, params):
""" use it for nearest interaction"""
d_12 = site_2.pos - site_1.pos
so = (1j*2/3)*params.lamda_rash*(sigma_x * d_12[1] - sigma_y * d_12[0])
return so
###
-------------------------------------------------------------------------
def Intrinsic_SOC(site_1, site_2, params):
""" use it for second nearest interaction"""
so = (+1j/(3*sqrt(3)))*params.lamda_Int*sigma_z
return so
###
-------------------------------------------------------------------------
def hopping(site_1, site_2, params):
return kinetic(site_1, site_2, params)+Rashba_SOC(site_1, site_2, params)
###
-------------------------------------------------------------------------
def hopping_lead(site_1, site_2, params):
return kinetic(site_1, site_2, params)+Rashba_SOC(site_1, site_2, params)
##
==========================================================================================================================
### ============ Make lattice ==============
params = SimpleNamespace()
params.lat_constant=0.24595
lat_constant = params.lat_constant
lat = kwant.lattice.honeycomb(a=lat_constant, norbs=2, name=["a", "b"])
a, b =lat.sublattices
def make_sys(l, w):
def rect(pos):
x, y = pos
return abs(x) < l and abs(y) < w
sys=kwant.Builder()
sys[lat.shape(rect, (1, 1))] = onsite
sys[lat.neighbors(1)] = hopping
sys.eradicate_dangling()
## clockwise sites for SOC (+1)
clocksites_a = (((-1, 0), a, a), ((0, 1), a, a), ((1, -1), a, a))
clocksites_b = (((1, 0), b, b), ((0, -1), b, b), ((-1, 1), b, b))
sys[[kwant.builder.HoppingKind(*hopping) for hopping in clocksites_a]]
= Intrinsic_SOC
sys[[kwant.builder.HoppingKind(*hopping) for hopping in clocksites_b]]
= Intrinsic_SOC
### creating leads
def lead_shape(pos):
x, y = pos
return abs(y)< w
sym = kwant.TranslationalSymmetry(lat.vec((-1, 0)))
sym.add_site_family(lat.sublattices[0], other_vectors=[(-1, 2)])
sym.add_site_family(lat.sublattices[1], other_vectors=[(-1, 2)])
lead = kwant.Builder(sym, conservation_law=-sigma_z)
lead[lat.shape(lead_shape, (1, 1))] = onsite
lead[lat.neighbors(1)] = hopping_lead
lead.eradicate_dangling()
# clockwise sites for SOC (+1)
clocksites_a = (((-1, 0), a, a), ((0, 1), a, a), ((1, -1), a, a))
clocksites_b = (((1, 0), b, b), ((0, -1), b, b), ((-1, 1), b, b))
lead[[kwant.builder.HoppingKind(*hopping) for hopping in clocksites_a]]
= Intrinsic_SOC
lead[[kwant.builder.HoppingKind(*hopping) for hopping in clocksites_b]]
= Intrinsic_SOC
sys.attach_lead(lead)
sys.attach_lead(lead.reversed())
return sys, lead
## ----------------- DIFINING TB AND SYSTEM PARAMETERS
----------------------
params = SimpleNamespace(Ep = 0, hop = -1, lat_constant= 0.24595,
center=(0, 0) )
## ----------------- PLOTTING THE SYSTEM ----------------------
syst, lead = make_sys(l=1, w=1)
fsyst= syst.finalized()
flead= lead.finalized()
# from ipywidgets import interact
# @interact(Int=(0, 1.5, .2), rash= (0, 2, .2), mass=(0, 0.51, .05),)
# def band_vs_Tranv(Int=0, rash =0, mass=0.,):
# fig = plt.figure(figsize=(7, 7))
# ax = fig.add_subplot(111)
# ax.set_xlabel("momentum [(lattice constant)^-1]", fontsize=15)
# ax.set_ylabel("Energy (eV) of the lead",fontsize=20)
# params.mass=mass
# params.lamda_Int = Int
# params.lamda_rash = rash
# kwant.plotter.bands(flead, show=False, ax =ax, args=[params])
params.mass=0.5
params.lamda_Int = 0.2
params.lamda_rash = 0.5
energy = params.mass+0.01
smatrix = kwant.smatrix(syst.finalized(),energy, args=[params])
total_trans = smatrix.transmission(1, 0)
up_up_trans = smatrix.transmission((1, 0), (0, 0))
print(total_trans)
print(up_up_trans)
dn_dn_trans = smatrix.transmission((1, 1), (0, 1))
up_dn_trans = smatrix.transmission((1, 1), (0, 0))
dn_up_trans = smatrix.transmission((1, 0), (0, 1))
print(total_trans)
print(up_up_trans)
print(dn_dn_trans)
print(up_dn_trans)
Le dim. 15 août 2021 à 15:30, Adel Belayadi <[email protected]> a écrit :
> Hey Anton,
> I should have thought of the physics behind my hopping values when calling
> the conservation law
> I have used realistic hopping values instead of np.eye(2, 6) and it is
> working now.
> Regards,
> Adel
>
> Le dim. 15 août 2021 à 14:43, Anton Akhmerov <[email protected]>
> a écrit :
>
>> Hi Adel,
>>
>> Running the code you provided also produces the following warning.
>>
>> > UserWarning: Hamiltonian breaks Conservation law, ignoring the symmetry
>> in the computation.
>>
>> This should explain what is happening.
>>
>> Best,
>> Anton
>>
>> On Sun, 15 Aug 2021 at 14:19, Adel Belayadi <[email protected]> wrote:
>> >
>> > Dear Anton;
>> >
>> > Thank you for your reply
>> >
>> > I have used the conservation_law as you have suggested. In fact, a
>> dictionary or function definition of the conservation law is working fine
>> now. I am not getting an error while defining the conservation law.
>> However, when calling smatrix.transmission((1, 1), (0, 1)) which stands
>> for spin up-up transmittance it raises: index 2 is out of bounds for axis 0
>> with size 2.
>> >
>> > Please see for example the following script which have two site
>> families, one with 2 orbitals and the second with 6 orbitals
>> >
>> > Cheers, Adel
>> >
>> >
>> > import kwant
>> > from matplotlib import pyplot
>> > import numpy as np
>> > from math import sqrt
>> > lat = kwant.lattice.honeycomb(a=1, norbs=[2, 6] , name=["a", "b"])
>> > a, b =lat.sublattices
>> > def make_sys(l, w, cons_law):
>> > def rect(pos):
>> > x, y = pos
>> > return abs(x) < l and abs(y) < w
>> > sys=kwant.Builder()
>> > sys[a.shape(rect, (0, 0))] = np.eye(2)
>> > sys[b.shape(rect, (0, 0))] = np.eye(6)
>> > ab_hopps= (((0, 0), a, b), ((0, 1), a, b), ((-1, 1), a, b))
>> > sys[[kwant.builder.HoppingKind(*hopping) for hopping in ab_hopps]]
>> = np.eye(2,6)
>> > sys.eradicate_dangling()
>> > sys[a.neighbors(1)] = np.eye(2)
>> > sys[b.neighbors(1)] = np.eye(6)
>> > ### creating leads
>> > def lead_shape(pos):
>> > x, y = pos
>> > return abs(y)< w
>> > sym = kwant.TranslationalSymmetry(lat.vec((-1, 0)))
>> > sym.add_site_family(lat.sublattices[0], other_vectors=[(-1, 2)])
>> > sym.add_site_family(lat.sublattices[1], other_vectors=[(-1, 2)])
>> > lead = kwant.Builder(sym, conservation_law=cons_law)
>> > lead[a.shape(rect, (0, 0))] = np.eye(2)
>> > lead[b.shape(rect, (0, 0))] = np.eye(6)
>> > lead[[kwant.builder.HoppingKind(*hopping) for hopping in ab_hopps]]
>> = np.eye(2,6)
>> > lead.eradicate_dangling()
>> > lead[a.neighbors(1)] = np.eye(2)
>> > lead[b.neighbors(1)] = np.eye(6)
>> > sys.attach_lead(lead)
>> > sys.attach_lead(lead.reversed())
>> > return sys
>> >
>> >
>> > syst = make_sys(l=3, w=2, cons_law= None)
>> > kwant.plot(syst ,fig_size=(8, 8))
>> > pyplot.show()
>> >
>> > print('************ calling the conservation law **************')
>> > sigma_z = np.array([[1, 0], [0, -1]])
>> > cons_law = {a: -sigma_z , b: -np.kron(sigma_z, np.eye(3))}
>> > syst = make_sys(l=3, w=2, cons_law=cons_law)
>> > smatrix = kwant.smatrix(syst.finalized(), 0.01)
>> > total_trans = smatrix.transmission(1, 0)
>> > up_up_trans = smatrix.transmission((1, 0), (0, 0))
>> > dn_dn_trans = smatrix.transmission((1, 1), (0, 1))
>> > up_dn_trans = smatrix.transmission((1, 1), (0, 0))
>> > dn_up_trans = smatrix.transmission((1, 0), (0, 1))
>> > print(total_trans)
>> > print(up_up_trans)
>> > print(dn_dn_trans)
>> > print(up_dn_trans)
>> > print(dn_up_trans)
>> > print(round(total_trans,
>> 3)==round(up_up_trans+dn_dn_trans+up_dn_trans+dn_up_trans, 3))
>> >
>> >
>> >
>> > Le sam. 14 août 2021 à 17:06, Anton Akhmerov <
>> [email protected]> a écrit :
>> >>
>> >> Hi Adel,
>> >>
>> >> Please check the documentation [1]. Since you have multiple site
>> >> families, you need to specify the conservation law as a dictionary, so
>> >> {metal_like_sub: -np.kron(sigma_z, np.eye(3)), C_like_sub1: -sigma_z,
>> >> C_like_sub2: -sigma_z}.
>> >>
>> >> Best,
>> >> Anton
>> >>
>> >>
>> https://kwant-project.org/doc/1/reference/generated/kwant.builder.Builder#kwant.builder.Builder
>> >>
>> >> On Sat, 14 Aug 2021 at 16:08, Adel Belayadi <[email protected]>
>> wrote:
>> >> >
>> >> > Dear Kwant user,
>> >> > Good day and I hope this email finds you well.
>> >> > Using conservation law was fully explained, in the Kwant mailing
>> list, while we are having sublattices with the same degree of freedom (same
>> number of orbitals). For instance, if we have two sites with 3 orbitals in
>> the case of a spin full Hamiltonian, the onsite matrix is
>> (onsite_matrix(6*6) for each site ). So in this case, the conservation law
>> will be straightforward: conservation_law=-np.kron(sigma_z, np.eye(3)).But
>> what would be the case if we are dealing with sublattices which have
>> different numbers of orbitals?
>> >> >
>> >> > In my case, I am interested in spin conductance. I am using 3
>> sublattices with different orbitals. Let us say for instance
>> >> > sublattice0 has (6-orbitals), sublattice1 with (2-orbitals) and
>> sublattice2 with (2-orbitals). In fact I am working on hetero-structure
>> with two sites of C-like(pz with spin up and down) and one site of
>> metal-like(3d-orbitals with spin up and down). I can resume my lattice as
>> follow
>> >> >
>> >> > lat = kwant.lattice.general(prim_vecs, basis=[ [...], [...], [...]]
>> , norbs=[6, 2, 2])
>> >> > metal-like_sub, C-like_sub1, C-like_sub2 = lat.sublattices
>> >> >
>> >> > According to my understanding the conservation law would be
>> something like
>> >> > conservation_law=-np.kron(sigma_z, np.eye(5)). However Kwant raises
>> an error with my choice of the conservation law.
>> >> >
>> >> > Is there something missing with my understanding of the conservation
>> law!. Or in the case of different orbitals per site, we need to define two
>> lattices (lattice1 for spin up and lattice2 for spindown).
>> >> > Thanks in advance,
>> >> > Best
>> >> >
>> >> >
>> >> >
>>
>