Lattice specification and bands

Checkerboard

Source code

"""Two dimensional checkerboard lattice with real hoppings"""
import pybinding as pb
import matplotlib.pyplot as plt
from math import pi

pb.pltutils.use_style()


def checkerboard(d=0.2, delta=1.1, t=0.6):
    lat = pb.Lattice(a1=[d, 0], a2=[0, d])
    lat.add_sublattices(
        ('A', [0, 0], -delta),
        ('B', [d/2, d/2], delta)
    )
    lat.add_hoppings(
        ([ 0,  0], 'A', 'B', t),
        ([ 0, -1], 'A', 'B', t),
        ([-1,  0], 'A', 'B', t),
        ([-1, -1], 'A', 'B', t)
    )
    return lat

lattice = checkerboard()
lattice.plot()
plt.show()
Checkerboard lattice, Brillouin zone and band structure
lattice.plot_brillouin_zone()
plt.show()
Checkerboard lattice, Brillouin zone and band structure
model = pb.Model(checkerboard(), pb.translational_symmetry())
solver = pb.solver.lapack(model)

bands = solver.calc_bands([0, 0], [0, 5*pi], [5*pi, 5*pi], [0, 0])
bands.plot()
plt.show()
Checkerboard lattice, Brillouin zone and band structure

Trestle

Source code

"""One dimensional lattice with complex hoppings"""
import pybinding as pb
import matplotlib.pyplot as plt

pb.pltutils.use_style()


def trestle(d=0.2, t1=0.8 + 0.6j, t2=2):
    lat = pb.Lattice(a1=1.3*d)
    lat.add_sublattices(
        ('A', [0,   0], 0),
        ('B', [d/2, d], 0)
    )
    lat.add_hoppings(
        (0, 'A', 'B', t1),
        (1, 'A', 'B', t1.conjugate()),
        (1, 'A', 'A', t2),
        (1, 'B', 'B', t2)
    )
    return lat


lattice = trestle()
lattice.plot()
plt.show()
Trestle lattice, Brillouin zone and band structure
lattice.plot_brillouin_zone()
plt.show()
Trestle lattice, Brillouin zone and band structure
model = pb.Model(trestle(), pb.translational_symmetry())
solver = pb.solver.lapack(model)

start, end = lattice.brillouin_zone()
bands = solver.calc_bands(start, end)
bands.plot()
plt.show()
Trestle lattice, Brillouin zone and band structure

Monolayer graphene

Source code

"""Create and plot a monolayer graphene lattice, its Brillouin zone and band structure"""
import pybinding as pb
import matplotlib.pyplot as plt
from math import sqrt, pi

pb.pltutils.use_style()


def monolayer_graphene():
    """Return the lattice specification for monolayer graphene"""
    a = 0.24595   # [nm] unit cell length
    a_cc = 0.142  # [nm] carbon-carbon distance
    t = -2.8      # [eV] nearest neighbour hopping

    # create a lattice with 2 primitive vectors
    lat = pb.Lattice(
        a1=[a, 0],
        a2=[a/2, a/2 * sqrt(3)]
    )

    lat.add_sublattices(
        # name and position
        ('A', [0, -a_cc/2]),
        ('B', [0,  a_cc/2])
    )

    lat.add_hoppings(
        # inside the main cell
        ([0,  0], 'A', 'B', t),
        # between neighboring cells
        ([1, -1], 'A', 'B', t),
        ([0, -1], 'A', 'B', t)
    )

    return lat


lattice = monolayer_graphene()
lattice.plot()
plt.show()
Monolayer graphene lattice, Brillouin zone and band structure
lattice.plot_brillouin_zone()
plt.show()
Monolayer graphene lattice, Brillouin zone and band structure
model = pb.Model(monolayer_graphene(), pb.translational_symmetry())
solver = pb.solver.lapack(model)

a_cc = 0.142
Gamma = [0, 0]
K1 = [-4*pi / (3*sqrt(3)*a_cc), 0]
M = [0, 2*pi / (3*a_cc)]
K2 = [2*pi / (3*sqrt(3)*a_cc), 2*pi / (3*a_cc)]

bands = solver.calc_bands(K1, Gamma, M, K2)
bands.plot(point_labels=['K', r'$\Gamma$', 'M', 'K'])
plt.show()
Monolayer graphene lattice, Brillouin zone and band structure
model.lattice.plot_brillouin_zone(decorate=False)
bands.plot_kpath(point_labels=['K', r'$\Gamma$', 'M', 'K'])
Monolayer graphene lattice, Brillouin zone and band structure

Monolayer graphene with Rashba spin-orbit coupling

Source code

"""Calculate the band structure of graphene with Rashba spin-orbit coupling"""
import pybinding as pb
import numpy as np
import matplotlib.pyplot as plt
from math import pi, sqrt


def monolayer_graphene_soc():
    """Return the lattice specification for monolayer graphene with Rashba SOC,
       see http://doi.org/10.1103/PhysRevB.95.165415 for reference"""
    from pybinding.constants import pauli
    from pybinding.repository.graphene import a_cc, a, t

    onsite = 0.05  # [eV] onsite energy
    rashba = 0.1   # [eV] strength of Rashba SOC
    rashba_so = 1j * 2/3 * rashba

    # create a lattice with 2 primitive vectors
    a1 = np.array([a / 2 * sqrt(3), a / 2])
    a2 = np.array([a / 2 * sqrt(3), -a / 2])
    lat = pb.Lattice(
        a1=a1, a2=a2
    )

    pos_a = np.array([-a_cc / 2, 0])
    pos_b = np.array([+a_cc / 2, 0])

    lat.add_sublattices(
        ('A', pos_a, [[ onsite, 0], [0,  onsite]]),
        ('B', pos_b, [[-onsite, 0], [0, -onsite]]))

    # nearest neighbor vectors
    d1 = (pos_b - pos_a) / a_cc          # [ 0,  0]
    d2 = (pos_b - pos_a - a1) / a_cc     # [-1,  0]
    d3 = (pos_b - pos_a - a2) / a_cc     # [ 0, -1]

    nn_hopp = np.array([[t, 0], [0, t]])                            # nn hopping, same spin
    t1 = nn_hopp + rashba_so * (pauli.x * d1[1] - pauli.y * d1[0])  # cross([sx , sy], [dx, dy])
    t2 = nn_hopp + rashba_so * (pauli.x * d2[1] - pauli.y * d2[0])
    t3 = nn_hopp + rashba_so * (pauli.x * d3[1] - pauli.y * d3[0])

    # name and position
    lat.add_hoppings(
        ([0,  0], 'A', 'B', t1),
        ([-1, 0], 'A', 'B', t2),
        ([0, -1], 'A', 'B', t3)
    )

    return lat


lattice = monolayer_graphene_soc()
lattice.plot()
plt.show()
Monolayer graphene with Rashba SOC, Brillouin zone and band structure
lattice.plot_brillouin_zone()
plt.show()
Monolayer graphene with Rashba SOC, Brillouin zone and band structure
model = pb.Model(lattice, pb.translational_symmetry())
solver = pb.solver.lapack(model)

k_points = model.lattice.brillouin_zone()
Gamma = [0, 0]
K1 = k_points[0]
K2 = k_points[2]
M = (k_points[0] + k_points[1]) / 2

bands = solver.calc_bands(K1, Gamma, M, K2)
bands.plot(point_labels=['K', r'$\Gamma$', 'M', 'K'])
plt.show()
Monolayer graphene with Rashba SOC, Brillouin zone and band structure

Monolayer graphene NN

Source code

"""Monolayer graphene with next-nearest hoppings"""
import pybinding as pb
import matplotlib.pyplot as plt
from math import sqrt, pi

pb.pltutils.use_style()


def monolayer_graphene_nn():
    a = 0.24595   # [nm] unit cell length
    a_cc = 0.142  # [nm] carbon-carbon distance
    t = -2.8      # [eV] nearest neighbour hopping
    t_nn = 0.25   # [eV] next-nearest neighbour hopping

    lat = pb.Lattice(
        a1=[a, 0],
        a2=[a/2, a/2 * sqrt(3)]
    )
    lat.add_sublattices(
        ('A', [0, -a_cc/2]),
        ('B', [0,  a_cc/2])
    )
    lat.add_hoppings(
        # between A and B inside the main cell
        ([0,  0], 'A', 'B', t),
        # between neighboring cells
        ([1, -1], 'A', 'B', t),
        ([0, -1], 'A', 'B', t),
        # next-nearest
        ([1,  0], 'A', 'A', t_nn),
        ([1,  0], 'B', 'B', t_nn),
        ([0,  1], 'A', 'A', t_nn),
        ([0,  1], 'B', 'B', t_nn),
        ([1, -1], 'A', 'A', t_nn),
        ([1, -1], 'B', 'B', t_nn)
    )
    return lat

lattice = monolayer_graphene_nn()
lattice.plot()
plt.show()
Monolayer graphene with next-nearest neighbor, Brillouin zone and band structure
lattice.plot_brillouin_zone()
plt.show()
Monolayer graphene with next-nearest neighbor, Brillouin zone and band structure
model = pb.Model(monolayer_graphene_nn(), pb.translational_symmetry())
solver = pb.solver.lapack(model)

a_cc = 0.142
Gamma = [0, 0]
K1 = [-4*pi / (3*sqrt(3)*a_cc), 0]
M = [0, 2*pi / (3*a_cc)]
K2 = [2*pi / (3*sqrt(3)*a_cc), 2*pi / (3*a_cc)]

# Note the elector-hole asymmetry in the band structure (due to t_nn).
bands = solver.calc_bands(K1, Gamma, M, K2)
bands.plot(point_labels=['K', r'$\Gamma$', 'M', 'K'])
plt.show()
Monolayer graphene with next-nearest neighbor, Brillouin zone and band structure

Bilayer graphene

Source code

"""Build the simplest model of bilayer graphene and compute its band structure"""
import pybinding as pb
import matplotlib.pyplot as plt
from math import sqrt, pi

pb.pltutils.use_style()


def bilayer_graphene():
    """Bilayer lattice in the AB-stacked form (Bernal-stacked)

    This is the simplest model with just a single intralayer and a single interlayer hopping.
    """
    a = 0.24595   # [nm] unit cell length
    a_cc = 0.142  # [nm] carbon-carbon distance
    c0 = 0.335    # [nm] interlayer spacing

    lat = pb.Lattice(a1=[a/2, a/2 * sqrt(3)], a2=[a/2, -a/2 * sqrt(3)])

    lat.add_sublattices(
        ('A1', [0,  -a_cc/2,   0]),
        ('B1', [0,   a_cc/2,   0]),
        ('A2', [0,   a_cc/2, -c0]),
        ('B2', [0, 3*a_cc/2, -c0])
    )

    lat.register_hopping_energies({
        'gamma0': -2.8,  # [eV] intralayer
        'gamma1': -0.4,  # [eV] interlayer
    })

    lat.add_hoppings(
        # layer 1
        ([ 0, 0], 'A1', 'B1', 'gamma0'),
        ([ 0, 1], 'A1', 'B1', 'gamma0'),
        ([-1, 0], 'A1', 'B1', 'gamma0'),
        # layer 2
        ([ 0, 0], 'A2', 'B2', 'gamma0'),
        ([ 0, 1], 'A2', 'B2', 'gamma0'),
        ([-1, 0], 'A2', 'B2', 'gamma0'),
        # interlayer
        ([ 0,  0], 'B1', 'A2', 'gamma1')
    )

    return lat


lattice = bilayer_graphene()
lattice.plot()
plt.show()
Simple bilayer graphene lattice, Brillouin zone and band structure
lattice.plot_brillouin_zone()
plt.show()
Simple bilayer graphene lattice, Brillouin zone and band structure
model = pb.Model(bilayer_graphene(), pb.translational_symmetry())
solver = pb.solver.lapack(model)

a_cc = 0.142
Gamma = [0, 0]
K1 = [-4*pi / (3*sqrt(3)*a_cc), 0]
M = [0, 2*pi / (3*a_cc)]
K2 = [2*pi / (3*sqrt(3)*a_cc), 2*pi / (3*a_cc)]

bands = solver.calc_bands(K1, Gamma, M, K2)
bands.plot(point_labels=['K', r'$\Gamma$', 'M', 'K'])
plt.show()
Simple bilayer graphene lattice, Brillouin zone and band structure

Phosphorene

Source code

"""Create and plot a phosphorene lattice, its Brillouin zone and band structure"""
import pybinding as pb
import matplotlib.pyplot as plt
from math import pi, sin, cos

pb.pltutils.use_style()


def phosphorene_4band():
    """Monolayer phosphorene lattice using the four-band model"""
    a = 0.222
    ax = 0.438
    ay = 0.332
    theta = 96.79 * (pi / 180)
    phi = 103.69 * (pi / 180)

    lat = pb.Lattice(a1=[ax, 0], a2=[0, ay])

    h = a * sin(phi - pi / 2)
    s = 0.5 * ax - a * cos(theta / 2)
    lat.add_sublattices(
        ('A', [-s/2,        -ay/2, h], 0),
        ('B', [ s/2,        -ay/2, 0], 0),
        ('C', [-s/2 + ax/2,     0, 0], 0),
        ('D', [ s/2 + ax/2,     0, h], 0)
    )

    lat.register_hopping_energies({
        't1': -1.22,
        't2': 3.665,
        't3': -0.205,
        't4': -0.105,
        't5': -0.055
    })

    lat.add_hoppings(
        # t1
        ([-1,  0], 'A', 'D', 't1'),
        ([-1, -1], 'A', 'D', 't1'),
        ([ 0,  0], 'B', 'C', 't1'),
        ([ 0, -1], 'B', 'C', 't1'),
        # t2
        ([ 0,  0], 'A', 'B', 't2'),
        ([ 0,  0], 'C', 'D', 't2'),
        # t3
        ([ 0,  0], 'A', 'D', 't3'),
        ([ 0, -1], 'A', 'D', 't3'),
        ([ 1,  1], 'C', 'B', 't3'),
        ([ 1,  0], 'C', 'B', 't3'),
        # t4
        ([ 0,  0], 'A', 'C', 't4'),
        ([ 0, -1], 'A', 'C', 't4'),
        ([-1,  0], 'A', 'C', 't4'),
        ([-1, -1], 'A', 'C', 't4'),
        ([ 0,  0], 'B', 'D', 't4'),
        ([ 0, -1], 'B', 'D', 't4'),
        ([-1,  0], 'B', 'D', 't4'),
        ([-1, -1], 'B', 'D', 't4'),
        # t5
        ([-1,  0], 'A', 'B', 't5'),
        ([-1,  0], 'C', 'D', 't5')
    )

    return lat

plt.figure(figsize=(6, 6))
lattice = phosphorene_4band()
lattice.plot()
plt.show()
Phosphorene lattice, Brillouin zone and band structure
lattice.plot_brillouin_zone()
plt.show()
Phosphorene lattice, Brillouin zone and band structure
model = pb.Model(phosphorene_4band(), pb.translational_symmetry())
solver = pb.solver.lapack(model)

ax = 0.438
ay = 0.332
kx = pi / ax
ky = pi / ay
bands = solver.calc_bands([kx, ky], [kx, 0], [0, 0], [0, ky], [kx, ky])
bands.plot(point_labels=["S", "Y", r"$\Gamma$", "X", "S"])
plt.show()
Phosphorene lattice, Brillouin zone and band structure
model.lattice.plot_brillouin_zone(decorate=False)
bands.plot_kpath(point_labels=["S", "Y", r"$\Gamma$", "X", "S"])
plt.show()
Phosphorene lattice, Brillouin zone and band structure