Dear all, thanks in advance for the support. As I already had occasion to explain to Xavier, my problem concerns a ring geometry where changing Hamiltonian values on nonexisting sites inside the ring (r<r1) actually changes the conductance behavior, while it should not. More precisely, I am building a ring, assigning to the ring sites a certain value given by the/onsite/ function, which however in its definition contains a conditional /if/instruction involving sites inside the ring. I was expecting that the /if/condition is then never satisfied and thus no change is made, but the conductance curve changes if I change the values of, say, the chemical potential (onsite energy) inside the ring, where no sites exist. However, If I do the /same/ changes /outside/ the ring (r>r2), no change occurs in the conductance, as it should. Of course one hypothesis is that I am in fact changing also Hamiltonian elements of the ring, but I don't see why. As suggested by Xavier, I will try to check that by using plotter. I attach the piece of code. I put a comment on the line where I write the conditional instruction for mu. If you try to run the code for different values of mu_inside, and/or different physical subregions inside the ring, you should see different conductance behavior. I am using python 2.7.5 on a Mac OS X 10.6.8 and the same holds true for a Mac OS X 10.8.5. Thanks a lot in advance for any hint on how to understand this misbehavior. Diego ============================================================================== import kwant import scipy, scipy.io from numpy import * from pylab import * import sys from numpy import mod from math import pi,atan2 import matplotlib as mpl from matplotlib import pyplot import random; import pickle params = {'backend': 'pdf', 'ytick.direction': 'out', 'text.usetex': True, 'font.size':16} rcParams.update(params) matplotlib.rcParams.update({'font.size': 20}) im=0+1j import tinyarray import scipy.sparse.linalg class SimpleNamespace(object): """A simple container for parameters.""" def __init__(self, **kwargs): self.__dict__.update(kwargs) s_0 = identity(2) s_z = array([[1, 0], [0, 1]]) s_x = array([[0, 1], [1, 0]]) s_y = array([[0, 1j], [1j, 0]]) #================================================================================= def make_system(a=1, W=10, r1=10, r2=20): lat = kwant.lattice.square(a) sys = kwant.Builder() def ring(pos): (x, y) = pos rsq = x ** 2 + y ** 2 return (r1 ** 2 < rsq < r2 ** 2) def mu(site, p): (x, y) = site.pos if abs(x) < r1 and abs(y) < r1: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring def onsite(site, p): return ( 4.0 * p.t  mu(site, p)) * s_0 def hopping(site0, site1, p): return p.t * s_0 # =============================================================================== sys[lat.shape(ring, (0, r1 + 1))] = onsite sys[lat.neighbors()] = hopping # =================================================================================== sym_lead = kwant.TranslationalSymmetry((a, 0)) lead = kwant.Builder(sym_lead) def lead_shape(pos): (x, y) = pos return (W / 2 < y < W / 2) lead[lat.shape(lead_shape, (0, 0))] = onsite lead[lat.neighbors()] = hopping sys.attach_lead(lead) sys.attach_lead(lead.reversed()) return sys def plot_conductance(sys, energies, p): data = [] for energy in energies: smatrix = kwant.smatrix(sys, energy, args=[p]) data.append(smatrix.transmission(1, 0)) axes([0.12, 0.12, 0.8, 0.8]) hold(True) pyplot.plot(energies, data, 'bo', ms=2.8, linewidth=0.5,markeredgewidth=0.03) pyplot.axis([0, 3.0, 0, 10.1]) pyplot.xlabel(r'$eV[t]$') pyplot.ylabel(r'$G[e^2/h]$') savefig('Ring_G(e).pdf') def main(): sys = make_system() kwant.plot(sys) sys = sys.finalized() params = SimpleNamespace(t=1.0, mu_ring=0.0, mu_inside=1.62) pyplot.figure() plot_conductance(sys, energies=[0.01 * i for i in xrange(300)], p=params) if __name__ == '__main__': main() ==================================================================================
Hi, Seems to me that in your onsite function:: def mu(site, p): (x, y) = site.pos if abs(x) < r1 and abs(y) < r1: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring the condition ``abs(x) < r1 and abs(y) < r1`` actually describes a *square* centred on (0, 0) with sides of length ``2 * r1``, and **not** a circle centred on (0, 0) with radius ``r1``. There are thus some sites at the corners which are included in this square, the onsites of which are affected when you change the value of ``mu_inside``. Changing the onsite function to the following should fix your problem:: def mu(site, p): (x, y) = site.pos if x**2 + y**2 < r1**2: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring Regards, Joe P.S. I have not tested the above code. You may need to change your shape function to:: def ring(pos): (x, y) = pos rsq = x ** 2 + y ** 2 return (r1 ** 2 <= rsq < r2 ** 2) # changed `<` to `<=` to avoid some possible annoying edge case
Thanks Joe for the comment. Actually that was just a stupid choice of subregion that I left at the end of my testing. The same problem is present if instead of the condition if abs(x) < r1 and abs(y) < r1: one uses if abs(x) < 2 and abs(y) < 2: or any similar small region /entirely/ included in the ring. Sorry for the useless extramistake that has misled already some people. Best, Diego On 24/3/14 11:47 PM, Joseph Weston wrote:
Hi,
Seems to me that in your onsite function::
def mu(site, p): (x, y) = site.pos if abs(x) < r1 and abs(y) < r1: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring
the condition ``abs(x) < r1 and abs(y) < r1`` actually describes a *square* centred on (0, 0) with sides of length ``2 * r1``, and **not** a circle centred on (0, 0) with radius ``r1``. There are thus some sites at the corners which are included in this square, the onsites of which are affected when you change the value of ``mu_inside``. Changing the onsite function to the following should fix your problem::
def mu(site, p): (x, y) = site.pos if x**2 + y**2 < r1**2: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring
Regards,
Joe
P.S. I have not tested the above code. You may need to change your shape function to::
def ring(pos): (x, y) = pos rsq = x ** 2 + y ** 2 return (r1 ** 2 <= rsq < r2 ** 2) # changed `<` to `<=`
to avoid some possible annoying edge case
Hi again, I've just had a thought as to what it could be, but I can't test it myself right now as I don't have access to a proper computer. I think this may be due to the way that kwant handles leads. When you create a Builder with a symmetry and add sites to it then kwant first translates the sites into the fundamental domain of the symmetry and adds those to the Builder. For translational symmetries the fundamental domain is the one which contains the point (0, 0). This means that when you construct your leads with:: lead[lat.shape(lead_shape, (0, 0))] = onsite then some of the onsites will actually evaluate to ``mu_inside`` as the lead sites satisfy the condition within ``onsite``. You can test this hypothesis by setting the "subregion" to something which does not include the lead's fundamental domain (in your case a single line of sites in the ydirection between W/2 and W/2 passing through the origin), or by setting the lead onsites "manually" for testing. Joe
Hi, thanks again for your time and your input. I tried what you said, by setting by hand the onsite elements of the lead, but it's not solving the problem. The result is identical to the one with the previous onsite function. Diego On 25/3/14 2:04 AM, Joseph Weston wrote:
Hi again,
I've just had a thought as to what it could be, but I can't test it myself right now as I don't have access to a proper computer. I think this may be due to the way that kwant handles leads.
When you create a Builder with a symmetry and add sites to it then kwant first translates the sites into the fundamental domain of the symmetry and adds those to the Builder. For translational symmetries the fundamental domain is the one which contains the point (0, 0). This means that when you construct your leads with::
lead[lat.shape(lead_shape, (0, 0))] = onsite
then some of the onsites will actually evaluate to ``mu_inside`` as the lead sites satisfy the condition within ``onsite``. You can test this hypothesis by setting the "subregion" to something which does not include the lead's fundamental domain (in your case a single line of sites in the ydirection between W/2 and W/2 passing through the origin), or by setting the lead onsites "manually" for testing.
Joe
Dear Diego, I tried yours script and I think that the issue is the one spotted by Joe: as leads are translationally invariant, you should also use functions to set their Hamiltonian matrix elements that are translationally invariant as Kwant does not guarantee the results otherwise (indeed internally the lead is translated to the fundamental domain but you don't really need to know that). The following script, adapted from yours, gives results which are independent from p.mu_inside if you use onsite_lead or dependent if you use onsite on the line where one sets the lead matrix elements Best regards, Xavier PS: the magic line to get the Hall resistance looks something like R= numpy.linalg.solve(G[:1,:1], [1, 0, 1]) where G is the transmission matrix in the Buttiker sense I = (e^2/h) G V PPS: Michael was faster than me, but I post this anyway as you can run the script below directly PPPS Joe, you're on vacation, remember? import kwant import scipy, scipy.io from numpy import * from pylab import * import sys from numpy import mod from math import pi,atan2 import matplotlib as mpl from matplotlib import pyplot import random; import pickle import tinyarray import scipy.sparse.linalg class SimpleNamespace(object): """A simple container for parameters.""" def __init__(self, **kwargs): self.__dict__.update(kwargs) s_0 = identity(2) s_z = array([[1, 0], [0, 1]]) s_x = array([[0, 1], [1, 0]]) s_y = array([[0, 1j], [1j, 0]]) #================================================================================= def make_system(a=1, W=10, r1=10, r2=20): lat = kwant.lattice.square(a) sys = kwant.Builder() def ring(pos): (x, y) = pos rsq = x ** 2 + y ** 2 return (r1 ** 2 < rsq < r2 ** 2) def mu(site, p): (x, y) = site.pos if abs(x) < r1/2 and abs(y) < r1/2: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring def onsite_lead(site, p): return ( 4.0 * p.t  p.mu_ring) * s_0 def onsite(site, p): return ( 4.0 * p.t  mu(site, p)) * s_0 def hopping(site0, site1, p): return p.t * s_0 # =============================================================================== sys[lat.shape(ring, (0, r1 + 1))] = onsite sys[lat.neighbors()] = hopping # =================================================================================== sym_lead = kwant.TranslationalSymmetry((a, 0)) lead = kwant.Builder(sym_lead) def lead_shape(pos): (x, y) = pos return (W / 2 < y < W / 2) lead[lat.shape(lead_shape, (0, 0))] = onsite # onsite_lead CHANGE THIS TO GET DIFFERENT BEHAVIOR lead[lat.neighbors()] = hopping sys.attach_lead(lead) sys.attach_lead(lead.reversed()) return sys def plot_conductance(sys, p): smatrix = kwant.smatrix(sys,3., args=[p]) print "CONDUCTANCE: ",smatrix.transmission(1, 0) def main(): sys = make_system() kwant.plot(sys) sys = sys.finalized() params = SimpleNamespace(t=1.0, mu_ring=0.0, mu_inside=1.02) plot_conductance(sys, p=params) if __name__ == '__main__': main() __________________________________________ Xavier Waintal SPSMSINACCEA,17 rue des Martyrs, 38054 Grenoble CEDEX 9, FRANCE Tel: 33 (0)4 38 78 03 27 email: xavier.waintal@cea.fr http://inac.cea.fr/Pisp/xavier.waintal __________________________________________ Le 25 mars 2014 à 02:38, diego <diego.rainis@unibas.ch> a écrit :
Hi,
thanks again for your time and your input.
I tried what you said, by setting by hand the onsite elements of the lead, but it's not solving the problem. The result is identical to the one with the previous onsite function.
Diego
On 25/3/14 2:04 AM, Joseph Weston wrote:
Hi again,
I've just had a thought as to what it could be, but I can't test it myself right now as I don't have access to a proper computer. I think this may be due to the way that kwant handles leads.
When you create a Builder with a symmetry and add sites to it then kwant first translates the sites into the fundamental domain of the symmetry and adds those to the Builder. For translational symmetries the fundamental domain is the one which contains the point (0, 0). This means that when you construct your leads with::
lead[lat.shape(lead_shape, (0, 0))] = onsite
then some of the onsites will actually evaluate to ``mu_inside`` as the lead sites satisfy the condition within ``onsite``. You can test this hypothesis by setting the "subregion" to something which does not include the lead's fundamental domain (in your case a single line of sites in the ydirection between W/2 and W/2 passing through the origin), or by setting the lead onsites "manually" for testing.
Joe
Dear all, thanks for all the answers. It is indeed due to the lead, I checked that now and it works if I separately define onsite for the lead without spacedependent conditions. That solves my problem, thanks a lot. Diego On 25/3/14 10:15 AM, Xavier Waintal wrote:
Dear Diego,
I tried yours script and I think that the issue is the one spotted by Joe: as leads are translationally invariant, you should also use functions to set their Hamiltonian matrix elements that are translationally invariant as Kwant does not guarantee the results otherwise (indeed internally the lead is translated to the fundamental domain but you don't really need to know that).
The following script, adapted from yours, gives results which are independent from p.mu_inside if you use onsite_lead or dependent if you use onsite on the line where one sets the lead matrix elements
Best regards,
Xavier
PS: the magic line to get the Hall resistance looks something like R= numpy.linalg.solve(G[:1,:1], [1, 0, 1]) where G is the transmission matrix in the Buttiker sense I = (e^2/h) G V
PPS: Michael was faster than me, but I post this anyway as you can run the script below directly
PPPS Joe, you're on vacation, remember?
import kwant import scipy, scipy.io from numpy import * from pylab import * import sys from numpy import mod from math import pi,atan2 import matplotlib as mpl from matplotlib import pyplot import random; import pickle
import tinyarray import scipy.sparse.linalg
class SimpleNamespace(object): """A simple container for parameters.""" def __init__(self, **kwargs): self.__dict__.update(kwargs) s_0 = identity(2) s_z = array([[1, 0], [0, 1]]) s_x = array([[0, 1], [1, 0]]) s_y = array([[0, 1j], [1j, 0]])
#=================================================================================
def make_system(a=1, W=10, r1=10, r2=20):
lat = kwant.lattice.square(a) sys = kwant.Builder()
def ring(pos): (x, y) = pos rsq = x ** 2 + y ** 2 return (r1 ** 2 < rsq < r2 ** 2)
def mu(site, p): (x, y) = site.pos if abs(x) < r1/2 and abs(y) < r1/2: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring
def onsite_lead(site, p): return ( 4.0 * p.t  p.mu_ring) * s_0 def onsite(site, p): return ( 4.0 * p.t  mu(site, p)) * s_0 def hopping(site0, site1, p): return p.t * s_0 # ===============================================================================
sys[lat.shape(ring, (0, r1 + 1))] = onsite sys[lat.neighbors()] = hopping # ===================================================================================
sym_lead = kwant.TranslationalSymmetry((a, 0)) lead = kwant.Builder(sym_lead) def lead_shape(pos): (x, y) = pos return (W / 2 < y < W / 2)
lead[lat.shape(lead_shape, (0, 0))] = onsite # onsite_lead CHANGE THIS TO GET DIFFERENT BEHAVIOR lead[lat.neighbors()] = hopping
sys.attach_lead(lead) sys.attach_lead(lead.reversed())
return sys
def plot_conductance(sys, p): smatrix = kwant.smatrix(sys,3., args=[p]) print "CONDUCTANCE: ",smatrix.transmission(1, 0)
def main(): sys = make_system() kwant.plot(sys) sys = sys.finalized() params = SimpleNamespace(t=1.0, mu_ring=0.0, mu_inside=1.02) plot_conductance(sys, p=params)
if __name__ == '__main__': main()
__________________________________________ Xavier Waintal SPSMSINACCEA,17 rue des Martyrs, 38054 Grenoble CEDEX 9, FRANCE Tel: 33 (0)4 38 78 03 27 email: xavier.waintal@cea.fr <mailto:xavier.waintal@cea.fr> http://inac.cea.fr/Pisp/xavier.waintal __________________________________________
Le 25 mars 2014 à 02:38, diego <diego.rainis@unibas.ch <mailto:diego.rainis@unibas.ch>> a écrit :
Hi,
thanks again for your time and your input.
I tried what you said, by setting by hand the onsite elements of the lead, but it's not solving the problem. The result is identical to the one with the previous onsite function.
Diego
On 25/3/14 2:04 AM, Joseph Weston wrote:
Hi again,
I've just had a thought as to what it could be, but I can't test it myself right now as I don't have access to a proper computer. I think this may be due to the way that kwant handles leads.
When you create a Builder with a symmetry and add sites to it then kwant first translates the sites into the fundamental domain of the symmetry and adds those to the Builder. For translational symmetries the fundamental domain is the one which contains the point (0, 0). This means that when you construct your leads with::
lead[lat.shape(lead_shape, (0, 0))] = onsite
then some of the onsites will actually evaluate to ``mu_inside`` as the lead sites satisfy the condition within ``onsite``. You can test this hypothesis by setting the "subregion" to something which does not include the lead's fundamental domain (in your case a single line of sites in the ydirection between W/2 and W/2 passing through the origin), or by setting the lead onsites "manually" for testing.
Joe
Elaborating on that: If one tries the same hamiltonian modification for different regions, one notices that: 1) for the inner square regions of the type if abs(x) < a and abs(y) < a: there is always an effect on the system, for any value of a, down to a=0. 2) Instead, surprisingly, if one makes modifications in the "outer" regions: if abs(x) > b and abs(y) > b: then there is no change in the system if b is BELOW r2, but above a certain threshold, in my system being b* = (r26). 3) Further, if instead one makes the changes in regions outside a circular domain, if x**2 + y**2 > (r22)**2: then THERE IS a modification to the system conductance. All this leads one to think that, maybe after the system finalization, the annulus sites have in fact been reshuffled and "packed" into a more convenient rectangular (square) shape, whose dimensions are determined by the total amount of sites. Here the annulus area is pi*(r2^2r1^2) ~ 950 sites, if one wants to pack them in a square, this gives a square side of 31, not far from my empirical estimate 2b*= 2(r26) = 28... (and the circle of radius (r22) tried above indeed crosses that square) Is it possible? On 24/3/14 11:47 PM, Joseph Weston wrote:
Hi,
Seems to me that in your onsite function::
def mu(site, p): (x, y) = site.pos if abs(x) < r1 and abs(y) < r1: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring
the condition ``abs(x) < r1 and abs(y) < r1`` actually describes a *square* centred on (0, 0) with sides of length ``2 * r1``, and **not** a circle centred on (0, 0) with radius ``r1``. There are thus some sites at the corners which are included in this square, the onsites of which are affected when you change the value of ``mu_inside``. Changing the onsite function to the following should fix your problem::
def mu(site, p): (x, y) = site.pos if x**2 + y**2 < r1**2: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring
Regards,
Joe
P.S. I have not tested the above code. You may need to change your shape function to::
def ring(pos): (x, y) = pos rsq = x ** 2 + y ** 2 return (r1 ** 2 <= rsq < r2 ** 2) # changed `<` to `<=`
to avoid some possible annoying edge case
Dear Diego, short answer is: No, kwant does not repack the sites after finalization  the coordinates of the sites stay where they were put in the beginning. To be honest, I don't fully understand now in which situations you see the dependence on p.mu_inside. But the problem can only be either in the system itself or in the leads. Let me give some examples on how to figure out where is the culprit: Check the onsite Hamiltonian of the system itself  In addition to plotting the system geometry, one can also plot any function defined on the system, and correspondingly also the onsite Hamiltonian elements. In your example, the easiest is to do this after finalization using kwant.plotter.map: vals = [sys.hamiltonian(i, i, params)[0, 0] for i in xrange(sys.graph.num_nodes)] kwant.plotter.map(sys, vals) (Note: this is for interacitve use, hence after removing the backend "pdf" from the beginning of your file) This plots the 0, 0 entry of your onsite matrix, and of course you can easily plot any entry you like. Using this plot you can see if you changed the onsite Hamiltonian somewhere in your system. As a side note, it is also possible to do this before finalization using kwant.plot: kwant.plot(sys, site_color=lambda site: sys[site](site, params)[0,0], site_symbol="s", site_size=0.5, cmap="jet") This just looks a bit more involved: with sys[site] you get the value at site, which in this case is a function (hence the function call afterwards) that returns a matrix (hence the [...]). Checking the leads  There is indeed one additional issue with your original code as noted by Joseph: You use the onsitefunction also in the lead. Now, the lead has translational invariance, and it is required (but not checked right now by kwant) that in this case also the onsitefunction must have translational invariance along the lead. The tricky thing here is that when we say translational invariance, this means not just going outwards in the lead: The lead is originally considered as infinite to both sides, and the Hamiltonian matrix elements must reflect this translational symmetry to both + and infinity. In practice this means that kwant takes some unit cell (with coordinates usually close to the origin of the coordinate system), and then repeats it to infinity. Now in your situation, this probably meant that the lead unit cell was affected by p.mu_inside. In your case this should be solved by simply assigning a constant value directly to the lead onsite values (instead of a function). However, you said this didn't solve your problem, so maybe there's something else. To check the leads, it is right now not easily possible to plot the onsite values as for the system. Instead, I advise you to check the resulting band structure of the leads using kwant.plotter.bands  if you see a difference in the band structures, you know you changed something in the leads. Best, Michael On 25.03.2014 00:35, diego wrote:
Elaborating on that:
If one tries the same hamiltonian modification for different regions, one notices that:
1) for the inner square regions of the type
if abs(x) < a and abs(y) < a:
there is always an effect on the system, for any value of a, down to a=0.
2) Instead, surprisingly, if one makes modifications in the "outer" regions:
if abs(x) > b and abs(y) > b:
then there is no change in the system if b is BELOW r2, but above a certain threshold, in my system being b* = (r26).
3) Further, if instead one makes the changes in regions outside a circular domain,
if x**2 + y**2 > (r22)**2:
then THERE IS a modification to the system conductance.
All this leads one to think that, maybe after the system finalization, the annulus sites have in fact been reshuffled and "packed" into a more convenient rectangular (square) shape, whose dimensions are determined by the total amount of sites. Here the annulus area is pi*(r2^2r1^2) ~ 950 sites, if one wants to pack them in a square, this gives a square side of 31, not far from my empirical estimate 2b*= 2(r26) = 28... (and the circle of radius (r22) tried above indeed crosses that square)
Is it possible?
On 24/3/14 11:47 PM, Joseph Weston wrote:
Hi,
Seems to me that in your onsite function::
def mu(site, p): (x, y) = site.pos if abs(x) < r1 and abs(y) < r1: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring
the condition ``abs(x) < r1 and abs(y) < r1`` actually describes a *square* centred on (0, 0) with sides of length ``2 * r1``, and **not** a circle centred on (0, 0) with radius ``r1``. There are thus some sites at the corners which are included in this square, the onsites of which are affected when you change the value of ``mu_inside``. Changing the onsite function to the following should fix your problem::
def mu(site, p): (x, y) = site.pos if x**2 + y**2 < r1**2: return p.mu_inside # this is the critical line: changing mu values inside the ring should not change the ring conductance else: return p.mu_ring
Regards,
Joe
P.S. I have not tested the above code. You may need to change your shape function to::
def ring(pos): (x, y) = pos rsq = x ** 2 + y ** 2 return (r1 ** 2 <= rsq < r2 ** 2) # changed `<` to `<=`
to avoid some possible annoying edge case
participants (5)

diego

Diego Rainis

Joseph Weston

Michael Wimmer

Xavier Waintal