Hello,
The kwant documentation on kwant.solvers.common.GreensFunction() states
that calling this method will return the retarded Green's function. Does
this method evaluate the equilibrium or non-equilibrium retarded Green's
function? I am interested in the time-dependent non-equilibrium density
matrix, which can be written as the kinetic lesser NE Green's function
G_<(t,t')
rho_neq (t) = -iG_<(t,t')|t=t'.
Thank you in advance.
Cheers,
Rudolf
Dear tkwant community,
Has anyone recently tried to simulate Josephson Junctions using tkwant?
I tried to follow the examples and tutorials, but several of them do not
work in the last tkwant version. I am sending attached and you can find
below the code that I am running but the current does not oscillate as a
function of time and the IV curve is linear. Any suggestions on what may be
wrong? Any help or comment is appreciated.
Here is the code:
import kwant
import tkwant
import tinyarray
import numpy as np
import matplotlib.pyplot as plt
import scipy.constants as ct
import time as comp_time
from mpi4py import MPI
comm=MPI.COMM_WORLD
start = comp_time.time()
tau_x = tinyarray.array([[0, 1], [1, 0]])
tau_y = tinyarray.array([[0, -1j], [1j, 0]])
tau_z = tinyarray.array([[1, 0], [0, -1]])
tau_0 = tinyarray.array([[1, 0], [0, 1]])
c_e = tinyarray.array([[1, 0], [0, 0]], complex)
c_h = tinyarray.array([[0, 0], [0, -1]], complex)
### parameters
Delta=0.02 #superconducting gap
barrier=0.9 #barrier in the N region
N=3 #size of the N region
S=10 #size of the S region
t=1
mu=t
I_time_avg=[]
V_values = np.linspace(-1.5*Delta, 1.5*Delta, 31)
def make_system(a=1):
lat = kwant.lattice.square(norbs=2)
syst = kwant.Builder()
#### Define junction. ####
### S
syst[(lat(x, 0) for x in range(-S, 0))] = (2 * t - mu) * tau_z + Delta
* tau_x
#### N
syst[(lat(x, 0) for x in range(0, N))] = (2 * t - mu + barrier) * tau_z
### S
syst[(lat(x, 0) for x in range(N, N+S))] = (2 * t - mu) * tau_z + Delta
* tau_x
# Hoppings
syst[lat.neighbors()] = -t * tau_z
### add JJ phase
syst[(lat(-1,0),lat(0,0))] = -t*np.exp(1j*phase(time))*c_e
-t*np.exp(-1j*phase(time))*c_h ### JJ phase
#### Define the leads. ####
lead_left = kwant.Builder(kwant.TranslationalSymmetry((-1, 0)))
lead_left[lat(0, 0)] = (2 * t - mu) * tau_z + Delta * tau_x
lead_left[lat.neighbors()] = -t * tau_z
lead_right = lead_left.reversed()
syst.attach_lead(lead_left)
syst.attach_lead(lead_right)
return syst
def phase(time):
### int_{0}^{t}A(1-cos(pi*t'/tau))dt'
def g(time,tau,A):
return A*time-A*tau/np.pi*np.sin(np.pi*time/tau)
if time<0:
return 0
elif time<tau:
return g(time, tau=tau, A=V/2)
else:
return V*(time-tau)
for V in V_values:
tau=0.01/(np.abs(0.04*np.pi*Delta))
tmax=4.01/(np.abs(0.04*np.pi*Delta))
syst = make_system()
lat = kwant.lattice.square(norbs=2)
# Finalize the system.
fsyst = syst.finalized()
#bondaries condition
boundaries=[tkwant.leads.MonomialAbsorbingBoundary(strength=8,degree=3,num_cells=1300)]
#kwant.plot(fsyst)
left_lead_interface = [(lat(1, 0), lat(0, 0))]
current_operator =
kwant.operator.Current(fsyst,where=left_lead_interface)
# occupation -- for each lead
occupation_l = tkwant.manybody.make_occupation(mu)
occupation_r = tkwant.manybody.make_occupation(mu+V)
occupation = [occupation_l, occupation_r]
# Create the solver
solver = tkwant.manybody.State(fsyst, tmax, occupation)
I_time = []
#get the current as a function of time
real_time = np.linspace(0, tmax, 200)
for time in real_time:
solver.evolve(time)
result = solver.evaluate(current_operator)
I_time.append(result)
data = np.column_stack((real_time, I_time))
I_time_avg.append(np.mean(I_time))
end = comp_time.time()
print('Execution time is: '+str((end - start)/60)+' minutes')
# check current as a function of time
plt.figure()
plt.plot(real_time, I_time)
#plot IV curve
plt.figure()
plt.plot(V_values, np.asarray(I_time_avg)/(2*barrier*Delta), 'o')
data = np.column_stack((V_values, I_time_avg))
np.savetxt('IV_curve.txt', data)
Best regards,
Denise
Hi Rudolf,
In tkwant the wavefunction approach is used instead of non-equilibrium Green functions. This provides better scalability in both time and system size. Non-equilibrium densities and currents can be obtained by evolving kwant operators in time. The Green functions (if needed) can be also obtained by integrating the wavefunctions over energy.
For more details you can have a look to the Weston’s thesis (book) [1] (mainly Chapter 2 and 3) as well as the tkwant tutorials.
Regards,
* J. Weston, Numerical Methods for Time-Resolved Quantum Nanoelectronics, Springer (2017)
Dear Joseph
At first thank you for kwant team. My problem for zigzag nanoribbon with
defect (hole) was solved with your help. Now I have a bit problem with
armchair nanoribbon. I have changed zigzag nanoribbon according to making
the armchair nanoribbon. It is OK when we have no defect on the system but
when I try to make hole on the scattering region and leads I have a
problem. The shape of holes on the scattering region and leads are not
same. My code is as the following. That is very king of you if you help me
where is the problem.
Thanks In advance for your help.
Beat wishes,
Nafise
import kwant
from math import sqrt
import matplotlib.pyplot as plt
import tinyarray
import numpy as np
import math
import cmath
#import scipy.linalg as la
import matplotlib
d=1.42;
a1=d*math.sqrt(3);
#on-site energy...................................................
t=-2.7;
latt = kwant.lattice.general([(0,a1),(a1*math.sqrt(3)/2,a1*0.5)],
[(-d/2,a1/2),(d/2,a1/2)])
a,b = latt.sublattices
syst= kwant.Builder()
def rectangle(pos):
x, y = pos
z=x**2+y**2
return -2.4*a1<x<2.4*a1 and -5*d<y<5*d
syst[latt.shape(rectangle,(1,1))]=0
def delet(pos):
x, y = pos
z=x**2+y**2
return z<(2*a1)**2
del syst[latt.shape(delet, (1,1))]
#nearest
neighbors.............................................................
syst[[kwant.builder.HoppingKind((0,0),a,b)]] =t
syst[[kwant.builder.HoppingKind((0,1),a,b)]] =t
syst[[kwant.builder.HoppingKind((-1,1),a,b)]] =t
ax=kwant.plot(syst);
sym = kwant.TranslationalSymmetry(latt.vec((-2,4)))
lead = kwant.Builder(sym)
def lead_shape(pos):
x, y = pos
return -5*d<y<5*d
lead[latt.shape(lead_shape, (1,1))]=0
def delet_lead(pos):
x, y = pos
z=x**2+y**2
return z<(2*a1)**2
del lead[latt.shape(delet_lead, (1,1))]
lead[[kwant.builder.HoppingKind((0,0),a,b)]]=t
lead[[kwant.builder.HoppingKind((0,1),a,b)]]=t
lead[[kwant.builder.HoppingKind((-1,1),a,b)]]=t
syst.attach_lead(lead,add_cells=0)
syst.attach_lead(lead.reversed(),add_cells=0)
ax=kwant.plot(syst);
Dear all,
I am using kwant to calculate the transmission of a bulk system. I use "
sys = kwant.wraparound.wraparound(sys) and lead =
kwant.wraparound.wraparound(lead, keep=0)" to change my tight-binding
ribbon system to a bulk system. The transmission is obtained using:
kys = np.arange(-0.5, 0.5, 0.001)
for ky in kys:
smatrix = kwant.smatrix(sys, energy, [ky])
transmission.append(smatrix.transmission(1, 0))
I want to obtain the Transmission vs the incident angle of an electron, and
I also need to do the integration ∫T(ky, E) dky or ∫T(θ,E)cosθdθ. I do not
undstand the ky in the system with kwant.wraparoun, what is the unit of ky,
what is k here for the system. Can you give me some advices on this?Thanks
in advance. The definition of incident angle is shown in the attached
Fig.1， the band structure for my bulk system is shown in attached Fig.2 and
the Tranmission vs ky is shown in Fig.3 (I do not know how the set the
values for ky).
Regards,
Hosein Khani
Dear Joseph
Thank you for your response. My problem for zigzag nanoribbon with defect
(hole) was solved with your help. Now I have a bit problem with armchair
nanoribbon. I have changed zigzag nanoribbon according to making the
armchair nanoribbon. It is OK when we have no defect on the system but when
I try to make hole on the scattering region and leads I have a problem. The
shape of holes on the scattering region and leads are not same. My code
is as the following. That is very king of you if you help me where is the
problem.
Thanks In advance for your help.
Beat wishes,
Nafise
import kwant
from math import sqrt
import matplotlib.pyplot as plt
import tinyarray
import numpy as np
import math
import cmath
#import scipy.linalg as la
import matplotlib
d=1.42;
a1=d*math.sqrt(3);
#on-site energy...................................................
t=-2.7;
latt = kwant.lattice.general([(0,a1),(a1*math.sqrt(3)/2,a1*0.5)],
[(-d/2,a1/2),(d/2,a1/2)])
a,b = latt.sublattices
syst= kwant.Builder()
def rectangle(pos):
x, y = pos
z=x**2+y**2
return -2.4*a1<x<2.4*a1 and -5*d<y<5*d
syst[latt.shape(rectangle,(1,1))]=0
def delet(pos):
x, y = pos
z=x**2+y**2
return z<(2*a1)**2
del syst[latt.shape(delet, (1,1))]
#nearest
neighbors.............................................................
syst[[kwant.builder.HoppingKind((0,0),a,b)]] =t
syst[[kwant.builder.HoppingKind((0,1),a,b)]] =t
syst[[kwant.builder.HoppingKind((-1,1),a,b)]] =t
ax=kwant.plot(syst);
sym = kwant.TranslationalSymmetry(latt.vec((-2,4)))
lead = kwant.Builder(sym)
def lead_shape(pos):
x, y = pos
return -5*d<y<5*d
lead[latt.shape(lead_shape, (1,1))]=0
def delet_lead(pos):
x, y = pos
z=x**2+y**2
return z<(2*a1)**2
del lead[latt.shape(delet_lead, (1,1))]
lead[[kwant.builder.HoppingKind((0,0),a,b)]]=t
lead[[kwant.builder.HoppingKind((0,1),a,b)]]=t
lead[[kwant.builder.HoppingKind((-1,1),a,b)]]=t
syst.attach_lead(lead,add_cells=0)
syst.attach_lead(lead.reversed(),add_cells=0)
ax=kwant.plot(syst);
Joseph Weston
[image: Attachments]Sun, Nov 3, 1:10 PM (1 day ago)
to me
Hi,
Please post Kwant usage questions to the Kwant mailing list, not my
personal email address.
Thanks,
Joe
On Sun, Nov 3, 2019 at 1:10 PM Joseph Weston <joseph.weston08(a)gmail.com>
wrote:
> Hi,
>
> Please post Kwant usage questions to the Kwant mailing list, not my
> personal email address.
>
> Thanks,
>
> Joe
>
> On 11/3/19 9:29 AM, Nafise Nouri wrote:
> > Dear Joseph
> >
> > At first thank you for your previous help. My problem for zigzag
> > nanoribbon with defect (hole) was solved with your help. Now I have a
> > bit problem with armchair nanoribbon. I have changed zigzag nanoribbon
> > according to making the armchair nanoribbon. It is OK when we have no
> > defect on the system but when I try to make hole on the scattering
> > region and leads I have a problem. The shape of holes on the
> > scattering region and leads are not same. My code is as the following.
> > That is very king of you if you help me where is the problem.
> >
> > Thanks In advance for your help.
> >
> > Beat wishes,
> > Nafise
> >
> > import kwant
> > from math import sqrt
> > import matplotlib.pyplot as plt
> > import tinyarray
> > import numpy as np
> > import math
> > import cmath
> > #import scipy.linalg as la
> > import matplotlib
> > d=1.42;
> > a1=d*math.sqrt(3);
> >
> > #on-site energy...................................................
> > t=-2.7;
> >
> > latt = kwant.lattice.general([(0,a1),(a1*math.sqrt(3)/2,a1*0.5)],
> > [(-d/2,a1/2),(d/2,a1/2)])
> > a,b = latt.sublattices
> > syst= kwant.Builder()
> >
> > def rectangle(pos):
> > x, y = pos
> > z=x**2+y**2
> > return -2.4*a1<x<2.4*a1 and -5*d<y<5*d
> >
> > syst[latt.shape(rectangle,(1,1))]=0
> >
> > def delet(pos):
> > x, y = pos
> > z=x**2+y**2
> > return z<(2*a1)**2
> >
> > del syst[latt.shape(delet, (1,1))]
> >
> > #nearest
> > neighbors.............................................................
> > syst[[kwant.builder.HoppingKind((0,0),a,b)]] =t
> > syst[[kwant.builder.HoppingKind((0,1),a,b)]] =t
> > syst[[kwant.builder.HoppingKind((-1,1),a,b)]] =t
> >
> > ax=kwant.plot(syst);
> >
> > sym = kwant.TranslationalSymmetry(latt.vec((-2,4)))
> >
> > lead = kwant.Builder(sym)
> >
> > def lead_shape(pos):
> > x, y = pos
> > return -5*d<y<5*d
> >
> > lead[latt.shape(lead_shape, (1,1))]=0
> >
> > def delet_lead(pos):
> > x, y = pos
> > z=x**2+y**2
> > return z<(2*a1)**2
> >
> > del lead[latt.shape(delet_lead, (1,1))]
> >
> > lead[[kwant.builder.HoppingKind((0,0),a,b)]]=t
> > lead[[kwant.builder.HoppingKind((0,1),a,b)]]=t
> > lead[[kwant.builder.HoppingKind((-1,1),a,b)]]=t
> >
> > syst.attach_lead(lead,add_cells=0)
> > syst.attach_lead(lead.reversed(),add_cells=0)
> > ax=kwant.plot(syst);
> >
> >
>
>
