raise RuntimeError("Numbers of left- and right-propagating " RuntimeError: Numbers of left- and right-propagating modes differ, possibly due to a numerical instability.
Hi all, I am getting the above error message while doing transport calculations in Silicene(next nearest neighbour hopping due to intrinsic SO). But it disappears when I set the the intrinsic SO to zero. What could be the possible reason. Thanks and Regards, Koustav Jana
Hi,
Could you please paste the part of the code that reproduces this error message?
Regards, Adel
On Sun, Jun 13, 2021 at 6:47 PM 170070051@iitb.ac.in wrote:
raise RuntimeError("Numbers of left- and right-propagating " RuntimeError: Numbers of left- and right-propagating modes differ, possibly due to a numerical instability.
Hi all, I am getting the above error message while doing transport calculations in Silicene(next nearest neighbour hopping due to intrinsic SO). But it disappears when I set the the intrinsic SO to zero. What could be the possible reason. Thanks and Regards, Koustav Jana
Hi, The code I am using is pasted below. This gives the above mentioned RunTime Error.
import kwant from math import pi,sqrt import numpy as np from matplotlib import pyplot from kwant.digest import uniform from types import SimpleNamespace
pauli_z = np.array([[1,0],[0,-1]]) sin_30, cos_30 = (1 / 2, sqrt(3) / 2) graphene = kwant.lattice.general([(1, 0), (sin_30, cos_30)], [(0, 0), (0, 1 / sqrt(3))])
a, b = graphene.sublattices
def lin0(y,W,jw) : if y < -jw : return -2 elif -jw <= y < jw : return 2*y/jw else : return 2
def lin1(y,W,jw) : if y < -W/6-jw : return -2 elif -W/6-jw <= y < -W/6+jw : return (y+W/6)/jw - 1 elif -W/6+jw <= y < W/6-jw : return 0 elif W/6-jw <= y < W/6+jw : return (y-W/6)/jw + 1 else : return 2
def lin2(y,W,jw) : if y < -jw : return 0 elif -jw <= y < jw : return y/jw+1 else : return 2
def make_system(W = 10*sqrt(3), L = 10, delta = 0, t = 1.6, lambda_so = 0) : lambda_so = 1j*lambda_so/(3*sqrt(3)) t_nn1_a = 1 * lambda_so * pauli_z # [ 1, 0] t_nn1_b = -1 * lambda_so * pauli_z t_nn2_a = -1 * lambda_so * pauli_z # [ 0, 1] t_nn2_b = 1 * lambda_so * pauli_z t_nn3_a = -1 * lambda_so * pauli_z # [ 1, -1] t_nn3_b = 1 * lambda_so * pauli_z
def channel(pos): x, y = pos return (0 <= x <= L) and (-W/2 < y <= W/2)
syst = kwant.Builder() del_fn = lambda y,W : lin1(y,W,W/20) def potential(site, U, U_disorder, Mex): (x, y) = site.pos salt = 0 d = -1 if (site.family == a) : d = 1 term1 = d*delta*del_fn(y,W)*np.eye(2) term2 = U*np.eye(2) term3 = Mex*pauli_z term4 = U_disorder * (uniform(repr(site), repr(salt)) - 0.5) * np.eye(2) return term1 + term2 + term3 + term4
def dummy(site, Mex): (x, y) = site.pos d = -1 if (site.family == a) : d = 1 term1 = d*delta*del_fn(y,W)*np.eye(2) term2 = Mex*pauli_z return term1 + term2
syst[graphene.shape(channel, (0, 0))] = potential hoppings = (((0, 0), a, b), ((0, 1), a, b), ((-1, 1), a, b)) syst[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = -t*np.eye(2) syst[kwant.builder.HoppingKind((1, 0), a, a)] = t_nn1_a syst[kwant.builder.HoppingKind((1, 0), b, b)] = t_nn1_b syst[kwant.builder.HoppingKind((0, 1), a, a)] = t_nn2_a syst[kwant.builder.HoppingKind((0, 1), b, b)] = t_nn2_b syst[kwant.builder.HoppingKind((1, -1), a, a)] = t_nn3_a syst[kwant.builder.HoppingKind((1, -1), b, b)] = t_nn3_b
# left lead sym0 = kwant.TranslationalSymmetry(graphene.vec((-1, 0)))
def lead0_shape(pos): x, y = pos return (-W/2 < y <= W/2)
lead0 = kwant.Builder(sym0) lead0[graphene.shape(lead0_shape, (0, 0))] = np.zeros((2,2)) lead0[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = -t*np.eye(2) lead0[kwant.builder.HoppingKind((1, 0), a, a)] = t_nn1_a lead0[kwant.builder.HoppingKind((1, 0), b, b)] = t_nn1_b lead0[kwant.builder.HoppingKind((0, 1), a, a)] = t_nn2_a lead0[kwant.builder.HoppingKind((0, 1), b, b)] = t_nn2_b lead0[kwant.builder.HoppingKind((1, -1), a, a)] = t_nn3_a lead0[kwant.builder.HoppingKind((1, -1), b, b)] = t_nn3_b
# right lead sym1 = kwant.TranslationalSymmetry(graphene.vec((1, 0)))
def lead1_shape(pos): x, y = pos return (-W/2 < y <= W/2)
lead1 = kwant.Builder(sym1) lead1[graphene.shape(lead1_shape, (0, 0))] = np.zeros((2,2)) lead1[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = -t*np.eye(2) lead1[kwant.builder.HoppingKind((1, 0), a, a)] = t_nn1_a lead1[kwant.builder.HoppingKind((1, 0), b, b)] = t_nn1_b lead1[kwant.builder.HoppingKind((0, 1), a, a)] = t_nn2_a lead1[kwant.builder.HoppingKind((0, 1), b, b)] = t_nn2_b lead1[kwant.builder.HoppingKind((1, -1), a, a)] = t_nn3_a lead1[kwant.builder.HoppingKind((1, -1), b, b)] = t_nn3_b
# dummy lead dum_sym = kwant.TranslationalSymmetry(graphene.vec((1, 0)))
def dum_lead_shape(pos): x, y = pos return (-W/2 < y <= W/2)
dum_lead = kwant.Builder(dum_sym) dum_lead[graphene.shape(dum_lead_shape, (0, 0))] = dummy dum_lead[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = -t*np.eye(2) dum_lead[kwant.builder.HoppingKind((1, 0), a, a)] = t_nn1_a dum_lead[kwant.builder.HoppingKind((1, 0), b, b)] = t_nn1_b dum_lead[kwant.builder.HoppingKind((0, 1), a, a)] = t_nn2_a dum_lead[kwant.builder.HoppingKind((0, 1), b, b)] = t_nn2_b dum_lead[kwant.builder.HoppingKind((1, -1), a, a)] = t_nn3_a dum_lead[kwant.builder.HoppingKind((1, -1), b, b)] = t_nn3_b
return syst, [lead0, lead1], dum_lead
W = 40*sqrt(3) L = 40 t = 1.3 E = 0.5 U = E lambda_so = 0.1 delta = -lambda_so U_disorder = 4*delta params = dict(U = U, U_disorder = U_disorder, Mex = 0)
syst, leads, dum_lead = make_system(W = W, L = L, delta = delta, t = t, lambda_so = lambda_so)
for lead in leads: syst.attach_lead(lead)
syst = syst.finalized()
local_dos = kwant.ldos(syst,energy=E,params=params) kwant.plotter.map(syst, local_dos[1::2], num_lead_cells=0, a=1/sqrt(3))
Hi, The following is the error message. This can help you find out where the problem is.
Traceback (most recent call last): File "/home/bmlabserver05/Koustav/kwant_simulations-kwant_valley_branch4/kwant_simulations-kwant_valley_branch1/kwant_valley_main/silicene_lattice.py", line 168, in <module> local_dos = kwant.ldos(syst,energy=E,params=params) File "/usr/lib/python3/dist-packages/kwant/solvers/common.py", line 535, in ldos params=params)[0] File "/usr/lib/python3/dist-packages/kwant/solvers/common.py", line 191, in _make_linear_sys prop, stab = lead.modes(energy, args, params=params) File "/usr/lib/python3/dist-packages/kwant/system.py", line 253, in modes return physics.modes(ham, hop, discrete_symmetry=symmetries) File "/usr/lib/python3/dist-packages/kwant/physics/leads.py", line 1143, in modes *symmetries) File "/usr/lib/python3/dist-packages/kwant/physics/leads.py", line 923, in compute_block_modes time_reversal, chiral) File "/usr/lib/python3/dist-packages/kwant/physics/leads.py", line 818, in make_proper_modes raise RuntimeError("Numbers of left- and right-propagating " RuntimeError: Numbers of left- and right-propagating modes differ, possibly due to a numerical instability.
Dear Sir or Madam, It seems that the number of modes from left and right are not equal while using the second nearest hopping. It should work for the first neighbor hopping but not further. You might overcome this issue by averting the use of a1 and a2 directions. You would simply add the following
sym0.add_site_family(graphene.sublattices[0], other_vectors=[(-1, 2)]) sym0.add_site_family(graphene.sublattices[1], other_vectors=[(-1, 2)])
sym1.add_site_family(graphene.sublattices[0], other_vectors=[(-1, 2)]) sym1.add_site_family(graphene.sublattices[1], other_vectors=[(-1, 2)])
dum_sym.add_site_family(graphene.sublattices[0], other_vectors=[(-1, 2)]) dum_sym.add_site_family(graphene.sublattices[1], other_vectors=[(-1, 2)])
I hope this will fix your issue. Best wishes Adel
Le lun. 14 juin 2021 à 05:37, 170070051@iitb.ac.in a écrit :
Hi, The following is the error message. This can help you find out where the problem is.
Traceback (most recent call last): File "/home/bmlabserver05/Koustav/kwant_simulations-kwant_valley_branch4/kwant_simulations-kwant_valley_branch1/kwant_valley_main/silicene_lattice.py", line 168, in <module> local_dos = kwant.ldos(syst,energy=E,params=params) File "/usr/lib/python3/dist-packages/kwant/solvers/common.py", line 535, in ldos params=params)[0] File "/usr/lib/python3/dist-packages/kwant/solvers/common.py", line 191, in _make_linear_sys prop, stab = lead.modes(energy, args, params=params) File "/usr/lib/python3/dist-packages/kwant/system.py", line 253, in modes return physics.modes(ham, hop, discrete_symmetry=symmetries) File "/usr/lib/python3/dist-packages/kwant/physics/leads.py", line 1143, in modes *symmetries) File "/usr/lib/python3/dist-packages/kwant/physics/leads.py", line 923, in compute_block_modes time_reversal, chiral) File "/usr/lib/python3/dist-packages/kwant/physics/leads.py", line 818, in make_proper_modes raise RuntimeError("Numbers of left- and right-propagating " RuntimeError: Numbers of left- and right-propagating modes differ, possibly due to a numerical instability.
Hi Adel, Thank you for your suggestion. This seems to work. One more request, could you please explain why this works. This seems to change the repeating unit in the leads, earlier it was inclined and now it has become vertical. What was wrong initially? Thanks
Dear Koustav,
Your example is puzzling me and I could not find the mistake! If you change the spin orbit coupling lambda_so = 0.1 to higher values like 0.6 and above it works! If you change the value of the width to 30 sqrt(1/3) or lower it works!
Please post the solution if you find it. Adel
On Mon, Jun 14, 2021 at 6:24 AM 170070051@iitb.ac.in wrote:
Hi, The code I am using is pasted below. This gives the above mentioned RunTime Error.
import kwant from math import pi,sqrt import numpy as np from matplotlib import pyplot from kwant.digest import uniform from types import SimpleNamespace
pauli_z = np.array([[1,0],[0,-1]]) sin_30, cos_30 = (1 / 2, sqrt(3) / 2) graphene = kwant.lattice.general([(1, 0), (sin_30, cos_30)], [(0, 0), (0, 1 / sqrt(3))])
a, b = graphene.sublattices
def lin0(y,W,jw) : if y < -jw : return -2 elif -jw <= y < jw : return 2*y/jw else : return 2
def lin1(y,W,jw) : if y < -W/6-jw : return -2 elif -W/6-jw <= y < -W/6+jw : return (y+W/6)/jw - 1 elif -W/6+jw <= y < W/6-jw : return 0 elif W/6-jw <= y < W/6+jw : return (y-W/6)/jw + 1 else : return 2
def lin2(y,W,jw) : if y < -jw : return 0 elif -jw <= y < jw : return y/jw+1 else : return 2
def make_system(W = 10*sqrt(3), L = 10, delta = 0, t = 1.6, lambda_so = 0) :
lambda_so = 1j*lambda_so/(3*sqrt(3)) t_nn1_a = 1 * lambda_so * pauli_z # [ 1, 0] t_nn1_b = -1 * lambda_so * pauli_z t_nn2_a = -1 * lambda_so * pauli_z # [ 0, 1] t_nn2_b = 1 * lambda_so * pauli_z t_nn3_a = -1 * lambda_so * pauli_z # [ 1, -1] t_nn3_b = 1 * lambda_so * pauli_z def channel(pos): x, y = pos return (0 <= x <= L) and (-W/2 < y <= W/2) syst = kwant.Builder() del_fn = lambda y,W : lin1(y,W,W/20) def potential(site, U, U_disorder, Mex): (x, y) = site.pos salt = 0 d = -1 if (site.family == a) : d = 1 term1 = d*delta*del_fn(y,W)*np.eye(2) term2 = U*np.eye(2) term3 = Mex*pauli_z term4 = U_disorder * (uniform(repr(site), repr(salt)) -
0.5) * np.eye(2) return term1 + term2 + term3 + term4
def dummy(site, Mex): (x, y) = site.pos d = -1 if (site.family == a) : d = 1 term1 = d*delta*del_fn(y,W)*np.eye(2) term2 = Mex*pauli_z return term1 + term2 syst[graphene.shape(channel, (0, 0))] = potential hoppings = (((0, 0), a, b), ((0, 1), a, b), ((-1, 1), a, b)) syst[[kwant.builder.HoppingKind(*hopping) for hopping in
hoppings]] = -t*np.eye(2) syst[kwant.builder.HoppingKind((1, 0), a, a)] = t_nn1_a syst[kwant.builder.HoppingKind((1, 0), b, b)] = t_nn1_b syst[kwant.builder.HoppingKind((0, 1), a, a)] = t_nn2_a syst[kwant.builder.HoppingKind((0, 1), b, b)] = t_nn2_b syst[kwant.builder.HoppingKind((1, -1), a, a)] = t_nn3_a syst[kwant.builder.HoppingKind((1, -1), b, b)] = t_nn3_b
# left lead sym0 = kwant.TranslationalSymmetry(graphene.vec((-1, 0))) def lead0_shape(pos): x, y = pos return (-W/2 < y <= W/2) lead0 = kwant.Builder(sym0) lead0[graphene.shape(lead0_shape, (0, 0))] = np.zeros((2,2)) lead0[[kwant.builder.HoppingKind(*hopping) for hopping in
hoppings]] = -t*np.eye(2) lead0[kwant.builder.HoppingKind((1, 0), a, a)] = t_nn1_a lead0[kwant.builder.HoppingKind((1, 0), b, b)] = t_nn1_b lead0[kwant.builder.HoppingKind((0, 1), a, a)] = t_nn2_a lead0[kwant.builder.HoppingKind((0, 1), b, b)] = t_nn2_b lead0[kwant.builder.HoppingKind((1, -1), a, a)] = t_nn3_a lead0[kwant.builder.HoppingKind((1, -1), b, b)] = t_nn3_b
# right lead sym1 = kwant.TranslationalSymmetry(graphene.vec((1, 0))) def lead1_shape(pos): x, y = pos return (-W/2 < y <= W/2) lead1 = kwant.Builder(sym1) lead1[graphene.shape(lead1_shape, (0, 0))] = np.zeros((2,2)) lead1[[kwant.builder.HoppingKind(*hopping) for hopping in
hoppings]] = -t*np.eye(2) lead1[kwant.builder.HoppingKind((1, 0), a, a)] = t_nn1_a lead1[kwant.builder.HoppingKind((1, 0), b, b)] = t_nn1_b lead1[kwant.builder.HoppingKind((0, 1), a, a)] = t_nn2_a lead1[kwant.builder.HoppingKind((0, 1), b, b)] = t_nn2_b lead1[kwant.builder.HoppingKind((1, -1), a, a)] = t_nn3_a lead1[kwant.builder.HoppingKind((1, -1), b, b)] = t_nn3_b
# dummy lead dum_sym = kwant.TranslationalSymmetry(graphene.vec((1, 0))) def dum_lead_shape(pos): x, y = pos return (-W/2 < y <= W/2) dum_lead = kwant.Builder(dum_sym) dum_lead[graphene.shape(dum_lead_shape, (0, 0))] = dummy dum_lead[[kwant.builder.HoppingKind(*hopping) for hopping in
hoppings]] = -t*np.eye(2) dum_lead[kwant.builder.HoppingKind((1, 0), a, a)] = t_nn1_a dum_lead[kwant.builder.HoppingKind((1, 0), b, b)] = t_nn1_b dum_lead[kwant.builder.HoppingKind((0, 1), a, a)] = t_nn2_a dum_lead[kwant.builder.HoppingKind((0, 1), b, b)] = t_nn2_b dum_lead[kwant.builder.HoppingKind((1, -1), a, a)] = t_nn3_a dum_lead[kwant.builder.HoppingKind((1, -1), b, b)] = t_nn3_b
return syst, [lead0, lead1], dum_lead
W = 40*sqrt(3) L = 40 t = 1.3 E = 0.5 U = E lambda_so = 0.1 delta = -lambda_so U_disorder = 4*delta params = dict(U = U, U_disorder = U_disorder, Mex = 0)
syst, leads, dum_lead = make_system(W = W, L = L, delta = delta, t = t, lambda_so = lambda_so)
for lead in leads: syst.attach_lead(lead)
syst = syst.finalized()
local_dos = kwant.ldos(syst,energy=E,params=params) kwant.plotter.map(syst, local_dos[1::2], num_lead_cells=0, a=1/sqrt(3))
Hi Abbout, Thanks for your reply. Even I had observed this, there was no error whenever I used narrow leads or a large value of the intrinsic SO. But the solution given by Adel Belayadi works for all the cases although I don't know why. You may have a look at it. Thanks