Thank you, Joseph! Your first advice was sufficient to make it 100 times faster. I strongly suspect that site.tag comparison was a very time-consuming operation. For the benefit of other users, I will add a subtlety that to arrange similar criterion-based hoppings for wires one also needs to compare with sites shifted by wire symmetry: e.g. sym = kwant.TranslationalSymmetry([1, 0]) def connectlead12(sys,crit,sym): ###only interlayer hoppings are included lat1 = (s for s in sys.sites() if s.family == lat) lat2 = (s for s in sys.sites() if s.family == lat2) pairs_of_sites = it.product(lat1, lat2) for site1, site2 in pairs_of_sites: if crit(site1, site2): yield (site1, site2) site2 = sym.act([1],site2) # act with the symmetry if crit(site1, site2): yield (site1, site2) Best wishes, Sergey On 25/05/16 15:12, Joseph Weston wrote:
Hi Sergey,
Dear Colleagues, Imagine, I have two different 2D lattices in a heterostructure and want to set up interlayer hoppings between the close-by atoms. I solve this problem in the way presented below. This works, but to my surprise it works so slowly that it became a bottle-neck of my whole computation. Is there a way to improve the code below? (I suspect, it may be the tag comparison that slows down the thing, but I am not sure)
def connect(sys,crit): for site1, site2 in it.product(sys.sites(), repeat=2): if crit(site1, site2): yield (site1, site2)
This loop will scale as O(N^2), so the first thing I would do is to try and cut down the number of sites N that you are comparing. For example, given that the sites are on different lattices, you could first create iterators that iterate over only the sites in each lattice independently:
lat_a_sites = (s for s in sys.sites() if s.family == lat_a) lat_b_sites = (s for s in sys.sites() if s.family == lat_b) pairs_of_sites = it.product(lat_a_sites, lat_b_sites): ...
This will already give you a fourfold speedup. If you can put in any other information that will allow you to reduce the number of sites that you compare (e.g. if your system is two 3D slabs, then only use sites "near" the interface)
I would then optimize the comparison in the following way:
def crit(site1, site2): return tinyarray.abs(site1.pos - site2.pos) < 0.1
sites_to_join = filter(pairs_of_sites, crit) sys[sites_to_join] = hopping
Let us know how this goes,
Happy Kwanting,
Joe