Dear John,
The solution I am proposing is not very clean, but at least it works fine.
To plot the wavefunction of a system for a fixed slice, you need to save
the result of the wavefunction for the desired sites, and then define a new
2D system recovering the slice you chose, and by reordering the results of
the wavefunction to fit the new order of the sites in the new system, the
plot becomes straightforward.
for your program it works like this:
# plots for the mode 0 and 1
sys = make_system_cyl(a=1, s=50, L = 15)
kwant.plot(sys)
sysf = sys.finalized()
params = SimpleNamespace(t=1)
wf = kwant.wave_function(sysf, 1.5, args=[params])
#define a new system for the desired slice.
lat2 = kwant.lattice.square()
sys_b = kwant.Builder() #new 2D system used for the plot
slice_x=5
mode=1
def Plot_mode(sys_b,slice_x,mode):
wavefunction = wf(1)[mode]
Sites= [site for (index,site) in enumerate(sysf.sites) if
site.pos[0]==slice_x ]
Indexes=[index for (index,site) in enumerate(sysf.sites) if
site.pos[0]==slice_x ]
waves=[wavefunction[index] for index in Indexes]
Sites2=[lat2(site.pos[1],site.pos[2]) for site in Sites] #used to
define the new system sys_b
sys_b[Sites2]=0
sys_b[lat2.neighbors()]=-1
Sites_b=[site for site in sys_b.finalized().sites] # sites of the new
system: with a different order!
Index_b=[Sites2.index(site) for site in Sites_b]
waves_b=[waves[index] for index in Index_b] # results is reordered
according to the new order of the sites.
return sys_b, waves_b
sys_b, waves_b= Plot_mode(sys_b=sys_b,slice_x=5,mode=0)
kwant.plotter.map(sys_b.finalized(),abs(array(waves_b))**2)
sys_b, waves_b= Plot_mode(sys_b=sys_b,slice_x=5,mode=1)
kwant.plotter.map(sys_b.finalized(),abs(array(waves_b))**2)
I hope this helps.
Regards,
Adel
On Tue, Apr 25, 2017 at 5:43 PM, John Eaves <johneave...@gmail.com> wrote:
> Good day,
>
> I'm trying to reconstruct something discussed in
> http://www.mail-archive.com/kwant-discuss@kwant-project.org/msg00432.html.
> In this answer it is described how to take the wavefunction from a 3D
> system and plot it in a 2D slice. This turned out to be something different
> than what the OP was asking for, but that is not so important. It is also
> mentioned that this has been asked before, but I couldn't find the right
> search terms. If my question has already been answered, I'd be very
> grateful for a referral to the relevant discussion!
>
> Now, as it turns out, I can't seem to get the above to work. I take some
> small cylinder, I calculate the wavefunction, but then I get stuck at how
> to extract for example the y = 0 plane and plot this in 2D. Below I'll post
> how far I got, perhaps the fix would be rather simple.
>
> def onsite(site, p):
> return 6*p.t
>
> def hopping(site0, site1, p):
> return -p.t
>
> def make_system_cyl(a=1, s=10, L = 15):
> lat = kwant.lattice.general([(0, 0, a), (a, 0, 0), (0, a, 0)])
>
>
> sys = kwant.Builder()
>
>
> def cyl_shape(pos):
> x, y, z = pos
> rsq = y ** 2 + z ** 2
> return rsq < s^2 and 0 <= x < L
>
>
> sys[lat.shape(cyl_shape, (1,1,1))] = onsite
> sys[lat.neighbors()] = hopping
>
>
> def lead_shape(pos):
> x, y, z = pos
> rsq = y ** 2 + z ** 2
> return rsq < s^2
>
> sym_lead = kwant.TranslationalSymmetry((-a,0,0))
> lead = kwant.Builder(sym_lead)
>
> lead[lat.shape(lead_shape, (-a,1,1))] = onsite
> lead[lat.neighbors()] = hopping
>
> sys.attach_lead(lead)
> sys.attach_lead(lead.reversed())
>
> return sys
>
> and I extract the wavefunction at some energy
>
> sys = make_system_cyl(a=1, s=10, L = 15)
> kwant.plot(sys)
> sysf = sys.finalized()
> params = SimpleNamespace(t=1)
> wf = kwant.wave_function(sysf, 1.5, args=[params])
>
> And then I try to extract the relevant parts
>
> wavefunction = wf(0)
> projected_wavefunction = np.empty((wavefunction.size, 1), dtype=np.complex)
>
> def in_sheet(index, site):
> return site.pos[1] == 0 # site in x-z plane
>
> sheet_idx_sites = filter(in_sheet, enumerate(sys.sites)
>
> for i, (idx, site) in enumerate(sheet_idx_sites):
> projected_wavefunction[i] = wavefunction[idx]
>
> Now this last part I took from the post referenced at the start, and it
> gives a syntax error. I'm not really sure how to fix it as it looks fine to
> me. What's wrong there?
>
> And more importantly, once this is fixed, how can one then make a
> kwant.plotter.map type of plot at this x-z plane?
>
--
Abbout Adel
