Posts How to Plot the First Brillouin Zone
Post
Cancel

How to Plot the First Brillouin Zone

Introduction

The first Brillouin Zone (BZ) is the Wigner-Seitz cell of the reciprocal lattice, which can be constructed by applying Voronoi decomposition to a lattice.


The reciprocal lattices (dots) and corresponding first Brillouin zones of (a) square lattice and (b) hexagonal lattice.

Constructing the Brillouin Zone

Below is the code that utilize scipy to generate the information of the BZ.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
#!/usr/bin/env python

import numpy as np

def get_brillouin_zone_3d(cell):
    """
    Generate the Brillouin Zone of a given cell. The BZ is the Wigner-Seitz cell
    of the reciprocal lattice, which can be constructed by Voronoi decomposition
    to the reciprocal lattice.  A Voronoi diagram is a subdivision of the space
    into the nearest neighborhoods of a given set of points. 

    https://en.wikipedia.org/wiki/Wigner%E2%80%93Seitz_cell
    https://docs.scipy.org/doc/scipy/reference/tutorial/spatial.html#voronoi-diagrams
    """

    cell = np.asarray(cell, dtype=float)
    assert cell.shape == (3, 3)

    px, py, pz = np.tensordot(cell, np.mgrid[-1:2, -1:2, -1:2], axes=[0, 0])
    points = np.c_[px.ravel(), py.ravel(), pz.ravel()]

    from scipy.spatial import Voronoi
    vor = Voronoi(points)

    bz_facets = []
    bz_ridges = []
    bz_vertices = []

    # for rid in vor.ridge_vertices:
    #     if( np.all(np.array(rid) >= 0) ):
    #         bz_ridges.append(vor.vertices[np.r_[rid, [rid[0]]]])
    #         bz_facets.append(vor.vertices[rid])

    for pid, rid in zip(vor.ridge_points, vor.ridge_vertices):
        # WHY 13 ????
        # The Voronoi ridges/facets are perpendicular to the lines drawn between the
        # input points. The 14th input point is [0, 0, 0].
        if(pid[0] == 13 or pid[1] == 13):
            bz_ridges.append(vor.vertices[np.r_[rid, [rid[0]]]])
            bz_facets.append(vor.vertices[rid])
            bz_vertices += rid

    bz_vertices = list(set(bz_vertices))

    return vor.vertices[bz_vertices], bz_ridges, bz_facets

The python function above returns the vertices, edges and facets of the BZ given the reciprocal basis vectors.

For an fcc lattice, the real-space basis vectors are defined by

1
2
3
cell = np.array([[0.0, 0.5, 0.5],
                 [0.5, 0.0, 0.5],
                 [0.5, 0.5, 0.0]])

For convenience, we have set the length of basis vectors to be 1.0. The reciprocal lattice basis vectors and lengths can then be obtained

1
2
3
4
5
# basis vectors of the reciprocal lattice
icell = np.linalg.inv(cell).T                

# length of the basis vectors
b1, b2, b3 = np.linalg.norm(icell, axis=1)   

The vertices, edges and facets of the BZ

1
v, e, f = get_brillouin_zone_3d(cell)

Plotting the Brillouin Zone

Once we have the information of the BZ, we can proceed to plot the Brillouin Zone either using Matplotlib or Mayavi.

Matplotlib

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

fig = plt.figure(
    figsize=(6, 6), dpi=300
)
ax = plt.subplot(111, projection='3d')
# ax = fig.add_subplot(111, projection='3d')

for xx in e:
    ax.plot(xx[:, 0], xx[:, 1], xx[:, 2], color='k', lw=1.0)

ax.set_xlim(-b1, b1)
ax.set_ylim(-b2, b2)
ax.set_zlim(-b3, b3)

plt.show()

FCC BZ Matplotlib

Mayavi

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
from mayavi import mlab
fig = mlab.figure(
    bgcolor=(1, 1, 1),
    size=(800, 800)
)

bz_line_width = b1 / 200

for xx in e:
    mlab.plot3d(xx[:, 0], xx[:, 1], xx[:, 2],
                tube_radius=bz_line_width,
                color=(0, 0, 0))

mlab.orientation_axes()
mlab.show()

FCC BZ Mayavi

The full code can be found here.

Plotly


FCC Brillouin Zone. Generated by Plotly.
This post is licensed under CC BY 4.0 by the author.