Dear Adel,
Thank you so much for your help.
I have some questions about it. Is the graphene "Finit" periodic in y
direction? Can I attach the leads to it?
I also have other questions related to current-voltage curve and applying
magnetic field. I wrote the following code for this aspect but I need to be
sure about it. Could you please help me?
Thank you in advance
###current-voltage####
def fermi_dirac(energy, mu, T):
return 1 / (np.exp((energy - mu) / (k_B * T)) + 1)
for voltage in voltages:
mu_left = voltage / 2
mu_right = -voltage / 2
# Compute the current using the Landauer-Büttiker formula
I = 0
# Loop over energy intervals and compute contributions
for e1, e2, g1, g2 in zip(energies[:-1], energies[1:],
conductance[:-1], conductance[1:]):
# Compute Fermi-Dirac distributions at e1 and e2 for left and
right reservoirs
f_left_e1 = fermi_dirac(e1, mu_left, temperature)
f_right_e1 = fermi_dirac(e1, mu_right, temperature)
f_left_e2 = fermi_dirac(e2, mu_left, temperature)
f_right_e2 = fermi_dirac(e2, mu_right, temperature)
# Current calculation using trapezoidal integration
delta_E = e2 - e1
I += ((g1 + g2) / 2) * ((f_left_e1 - f_right_e1) + (f_left_e2 -
f_right_e2)) / 2 * delta_E
# Convert from eV to Amperes
current.append(I * e/h * 1e6) # Multiply by e to convert the
current to Amperes
return current
###magnetic field###
def hopping_B(site1, site2, phi):
x1,y1=site1.pos
x2,y2=site2.pos
y = 0.5 * (y1 + y2)
y *= 1e-10
A_x = phi * y
peierls = A_x * (x1 - x2)
peierls *= 1e-10
return -t * np.exp(1j * 2*np.pi/phi0 * peierls)
On Thu, 14 Nov 2024, 10:43 Adel Belayadi, <[email protected]> wrote:
> Dear Masoumeh,
> I share an example of graphene with 04 atoms in the unit cell. You can
> also have the zigzag and armchair ribbon from that.
> I hope this helps
> Best
> Adel
>
> import kwant
> import matplotlib.pyplot as plt
> from math import sqrt
>
> def graphene_04_atoms (a=1.42, t=-1, Ep=0, Lx=5, Ly=5, boundary ="zigzag"):
> a1 = a*sqrt(3)
> a2 = 3*a
>
> ### here premitive vectors are perpondicular:
> primitive_vector_1.primitive_vector_1 = 0
> primitive_vector_1 = (a1, 0)
> primitive_vector_2 = (0, a2)
>
> ## Now let us define the atoms positions of the unit cell
> Pos_A1 = ( 0, -a/2)
> Pos_B1 = ( 0, a/2)
> Pos_A2 = ( a1/2, a)
> Pos_B2 = ( a1/2, 2*a)
>
>
> lat = kwant.lattice.general([primitive_vector_1, primitive_vector_2],
> [Pos_A1, Pos_B1, Pos_A2,
> Pos_B2])
> Sub_A1, Sub_B1, Sub_A2, Sub_B2 = lat.sublattices
>
>
> #...................................................................................
> if boundary == "Finit":
> syst= kwant.Builder()
> def square(pos):
> x, y = pos
> return abs(x)<Lx and abs(y)<Ly
> syst[lat.shape(square, (0,0))]= Ep
> #nearest
> neighbors.............................................................
> ## Be carfull if you consider second nearest neighbors
> ## =========================================================
> ## Hopping within unit cell ================================
> syst[[kwant.builder.HoppingKind((0,0),Sub_A1,Sub_B1)]] = t
> syst[[kwant.builder.HoppingKind((0,0),Sub_B1,Sub_A2)]] = t
> syst[[kwant.builder.HoppingKind((0,0),Sub_A2,Sub_B2)]] = t
> ## Hopping between neighbouring unit cell=======================
> syst[[kwant.builder.HoppingKind((+1, +1),Sub_A1,Sub_B2)]] = t
> syst[[kwant.builder.HoppingKind(( 0, +1),Sub_A1,Sub_B2)]] = t
> syst[[kwant.builder.HoppingKind((+1, 0),Sub_B1,Sub_A2)]] = t
>
> syst.eradicate_dangling()
> # Plot the closed system without leads.
> kwant.plot(syst);
>
> if boundary == "zigzag":
> def square(pos):
> x, y = pos
> return abs(y)<Ly
> sym = kwant.TranslationalSymmetry(primitive_vector_1)
> syst = kwant.Builder(sym)
> syst[lat.shape(square, (0,0))]= Ep
> syst[[kwant.builder.HoppingKind((0,0),Sub_A1,Sub_B1)]] = t
> syst[[kwant.builder.HoppingKind((0,0),Sub_B1,Sub_A2)]] = t
> syst[[kwant.builder.HoppingKind((0,0),Sub_A2,Sub_B2)]] = t
> ## Hopping between neighbouring unit cell=======================
> syst[[kwant.builder.HoppingKind((+1, +1),Sub_A1,Sub_B2)]] = t
> syst[[kwant.builder.HoppingKind(( 0, +1),Sub_A1,Sub_B2)]] = t
> syst[[kwant.builder.HoppingKind((+1, 0),Sub_B1,Sub_A2)]] = t
>
> syst.eradicate_dangling()
> # Plot the closed system without leads.
> kwant.plot(syst);
>
> if boundary == "armchair":
> def square(pos):
> x, y = pos
> return abs(x)<Lx
> sym = kwant.TranslationalSymmetry(primitive_vector_2)
> syst = kwant.Builder(sym)
> syst[lat.shape(square, (0,0))]= Ep
> syst[[kwant.builder.HoppingKind((0,0),Sub_A1,Sub_B1)]] = t
> syst[[kwant.builder.HoppingKind((0,0),Sub_B1,Sub_A2)]] = t
> syst[[kwant.builder.HoppingKind((0,0),Sub_A2,Sub_B2)]] = t
> ## Hopping between neighbouring unit cell=======================
> syst[[kwant.builder.HoppingKind((+1, +1),Sub_A1,Sub_B2)]] = t
> syst[[kwant.builder.HoppingKind(( 0, +1),Sub_A1,Sub_B2)]] = t
> syst[[kwant.builder.HoppingKind((+1, 0),Sub_B1,Sub_A2)]] = t
>
> syst.eradicate_dangling()
> # Plot the closed system without leads.
> kwant.plot(syst);
>
>
> graphene_04_atoms (a=1.42, t=-1, Ep=0, Lx=5, Ly=15, boundary ="Finit");
> graphene_04_atoms (a=1.42, t=-1, Ep=0, Lx=5, Ly=15, boundary ="zigzag");
> graphene_04_atoms (a=1.42, t=-1, Ep=0, Lx=15, Ly=5, boundary ="armchair");
>
>
> Le lun. 11 nov. 2024 à 14:20, Masoumeh Alihosseini via Kwant-discuss <
> [email protected]> a écrit :
>
>> Hi dear all,
>> I have a question, I want to work with rectangular graphene cell with 4
>> atoms in the unitcell. Like zigzag nanoribbon but I want to apply
>> periodicity in y direction for both scattering region and leads. Could you
>> please help me?
>> Thank you,
>> Masoumeh
>>
>