To add to the discussion, Joe Weston recently implemented a cleaned up version of the bound state solver, available here: https://gitlab.kwant-project.org/kwant/kwant/merge_requests/320
If you don't want to install the development version of Kwant, just grab the boundstate.py file from that merge request. Best, Anton On Tue, 10 Dec 2019 at 09:38, Xavier Waintal <xavier.wain...@cea.fr> wrote: > > Dear Denise, > > Slowly decaying bound states can happen whenever the bound state energy is > close to the opening > of a channel (e.g. close to the gap in your case). To move it away from the > gap, you can put a phase difference across the superconductors. > > As a sanity check you can take a finite system and perform an exact > diagonalisation (for several sizes that include a finite portion of the > leads). This is imprecise when the bound states are slowly decaying, but it > should allow you to spot obvious problems. > > Note that the algorithm is still experimental. As it requires some root > finding for a function f(E)=0, it is always a good idea to plot the function > f(E) to see if the algorithm has missed a root. > > Hope this helps. > > Best regards, > Xavier > > > Le 9 déc. 2019 à 19:27, Denise Puglia <dpuglia....@gmail.com> a écrit : > > Dear All, > > Has anyone tried to calculate bound states in a JJ? I used the package by > Benoit Rossignol available at: > https://gitlab.kwant-project.org/kwant/boundstate . > I am sending attached the code I used to simulate it but I am having some > difficulties interpreting the results. First, the algorithm could not find > bound states, i.e., psi_tot returned by the solver is an empty array, for > N<5. Second, it returns a total wavefunction for larger N, however the > wavefunction oscillates between each site and decays very slowly in the > leads. In the original proposal of this method, Istas et al > (arXiv:1711.08250) solves the wavefunction of quantum billiard (Fig 5) in > which this does not seem to happen. Has anyone found something similar? I do > know if it a mistake in the code, an attribute of discretization or something > else. > > Scanning over N=[1,2,3,5,10,20] and ,u=[2, Delta, 0, -2] the algorithm > returns the following non-zero wavefunctions: > SNS junction with N=10, S=20 and mu=2 > <N_10_S_20_mu_2.png> > SNS junction with N=20, S=20 and mu=2 > <N_20_S_20_mu_2.png> > SNS junction with N=5, S=20 and mu=-2 (in this case the amplitude of the > bound states decays to zero in the lead) > <N_5_S_20_mu_-2.png> > > SNSNS junction with 20/10/10/10/20 sites and mu=2. The wavefunction in the N > region is similar to the SNS junction but there seems to be no decay in the > middle S region. > <N_10_S_20_mu_2.png> > and of course it's mirrored version: > <NSN_10_mu_2.png> > > I appreciate any comments on this subject. > > Best regards, > Denise > > <BS_calc.py><fast.py><_common.py> > >
