Re: [Kwant] Wavefunction for 2D NISIN

[Joseph Weston]

Hi,


>  
> I still don’t understand exactly how I can implement this. From what I
> can understand you choose the energy not the eigenstate for the
> density function.


I am unclear what you mean here. In your code you do the following:


    ham_mat = sys.hamiltonian_submatrix(sparse=True, params = params)
    evals, evecs = sorted_eigs(sla.eigsh(ham_mat.tocsc(), k=20, sigma=0))

    kwant.plotter.map(sys, np.abs(evecs[:, 9])**2


You take the hamiltonian of the scattering region, diagonalize it, and
then plot *a single eigenvector*. Your problem was that because you have
>1 degree of freedom per site you cannot simply plot the eigenvector
as-is. You need to first compute the density per site for this
eigenvector, which you can do using the operator module:


    rho = kwant.operator.Density(np.kron(-s_z, s_0))  # holes count +1
to the density, and electrons -1

    kwant.plotter.map(sys, rho(evecs[:, 9))


Note that you call the density operator with a single eigenvector.


Independently of the above it's not clear to me that you should be
diagonalizing just the scattering region Hamiltonian, given that your
system has leads. It seems to me that you should be calculating the
scattering states and the density of those.


Happy Kwanting,


Joe




signature.asc
Description: OpenPGP digital signature

Re: [Kwant] Wavefunction for 2D NISIN

[Tidjani, Hanifa]

Hi Anton,

I still don’t understand exactly how I can implement this. From what I can 
understand you choose the energy not the eigenstate for the density function.

Is there a way to use the eigenvector solution in 2.4 for the wavefunction for 
a system with two sublattices?

Or am I simply misunderstanding where in the density function you encode the 
eigenstate?



Kind regards,

Hanifa

From: Anton Akhmerov 
Sent: Thursday, August 29, 2019 2:04:47 PM
To: Tidjani, Hanifa 
Cc: kwant-discuss@kwant-project.org 
Subject: Re: [Kwant] Wavefunction for 2D NISIN

Hi Hanifa,

This is certainly possible. It's up to you for which wave function you
compute the density. You can compute the eigenstates and compute the
expectation all the same.

Best,
Anton

On Mon, 26 Aug 2019 at 18:46, Tidjani, Hanifa  wrote:
>
> Dear Anton,
>
>
>
> Thank you for your reply.
>
>
>
> This does indeed work, however I am concerned that I am driving it at a 
> certain energy and then plotting the wave function. Is there a way to plot 
> the wave function of a particular eigenstate rather than at a particular 
> energy as is done in tutorial 2.4?
>
>
>
> Kind regards,
>
>
>
> Hanifa
>
>
>
> 
> From: Anton Akhmerov 
> Sent: Monday, August 26, 2019 9:16:59 AM
> To: Tidjani, Hanifa 
> Cc: kwant-discuss@kwant-project.org 
> Subject: Re: [Kwant] Wavefunction for 2D NISIN
>
> Hi Hanifa,
>
> The mismatch between the eigenvector size and the number of sites in
> the system is because your system has more than one degree of freedom
> per site. Please see how to deal with such cases in the tutorial:
> https://eur03.safelinks.protection.outlook.com/?url=https%3A%2F%2Fkwant-project.org%2Fdoc%2F1%2Ftutorial%2Foperatorsdata=01%7C01%7Chanifa.tidjani%40kcl.ac.uk%7Cfe0e226f1cb240cc912608d72c817fe2%7C8370cf1416f34c16b83c724071654356%7C0sdata=Nr9oHbH%2FRzVmWJaYBwpYapAaXO3I%2BZZTIR6ui3yWw0g%3Dreserved=0
>
> Best,
> Anton
>
> On Mon, 26 Aug 2019 at 10:08, Tidjani, Hanifa  
> wrote:
> >
> >
> >
> > I am trying to plot the wavefunction for a 2D NISIN Hamiltonian however I 
> > am unable to get a result, this is my code:
> >
> >
> >
> >
> >
> > import kwant
> >
> > import tinyarray as ta
> >
> > import numpy as np
> >
> > import matplotlib.pyplot as plt
> >
> > from tqdm import tqdm
> >
> > import winsound
> >
> > import os
> >
> > from datetime import datetime
> >
> > import scipy.sparse.linalg as sla
> >
> >
> >
> >
> >
> > sig_0 = np.identity(4)
> >
> > s_0 = np.identity(2)
> >
> > s_z = np.array([[1., 0.], [0., -1.]])
> >
> > s_x = np.array([[0., 1.], [1., 0.]])
> >
> > s_y = np.array([[0., -1j], [1j, 0.]])
> >
> >
> >
> > tau_z = ta.array(np.kron(s_z, s_0))
> >
> > tau_x = ta.array(np.kron(s_x, s_0))
> >
> > tau_y = ta.array(np.kron(s_y, s_0))
> >
> > sig_z = ta.array(np.kron(s_0, s_z))
> >
> > tau_zsig_x = ta.array(np.kron(s_z, s_x))
> >
> > tau_zsig_y = ta.array(np.kron(s_z, s_y))
> >
> >
> >
> > # Lorentzian definition
> >
> >
> >
> >
> >
> > def lorentzianxy(x, y, ind,y0):
> >
> > gam = 3
> >
> > L = params["L"]
> >
> > #if 0 >
> > return gam ** 2 /( gam ** 2 + (x - ind) ** 2 + (y - y0) ** 2)
> >
> ># else :
> >
> > #return 0
> >
> >
> >
> > # Define potential
> >
> >
> >
> > def potential(site, pot,ind,y0):
> >
> > (x, y) = site.pos
> >
> > return pot * lorentzianxy(x,y,ind,y0)
> >
> >
> >
> >
> >
> > # Make system
> >
> >
> >
> > def makeNISIN2D(params):
> >
> >
> >
> > # List all the params
> >
> > W=params["W"]
> >
> > L=params["L"]
> >
> > pot = params["pot"]
> >
> > ind = params["ind"]
> >
> > mu=params["mu"]
> >
> > B=params["B"]
> >
> > Delta=params["Delta"]
> >
> > alpha=params["alpha"]
> >
> > t=params["t"]
> >
> > barrier = params["barrier"]
> >
> > Smu = params["Smu"]
> >
> >
> >
> >
> >
> >
> >
> > def lorentzianxy(x, y, ind,y0):
> >
> > gam = 3
> >
> > L = params["L"]
> >
> > #if 0 >
> > return gam ** 2 /( gam ** 2 + (x - ind) ** 2 + (y - y0) ** 2)
> >
> ># else :
> >
> > #return 0
> >
> >
> >
> > # Define potential
> >
> >
> >
> > def potential(site, pot,ind,y0):
> >
> > (x, y) = site.pos
> >
> > return pot * lorentzianxy(x,y,ind,y0)
> >
> > #if y < lorentzian(x,ind)*4:
> >
> > #
> >
> > #return (tru_pot)
> >
> > #else:
> >
> > #return 0
> >
> > #
> >
> > def Straightpotential(site, pot,ind):
> >
> > (x, y) = site.pos
> >
> > if x == ind:
> >
> > return pot
> >
> > else:
> >
> > return 0
> >
> >
> >
> > def onsiteSc(site, pot,ind,y0):
> >
> >