Quasi-crystal

In this tutorial, we show how to construct the quasi-crystal, in which we also need to shift twist, reshape and merge the cells. Taking bilayer graphene quasicrystal as an example, a quasicrystal with 12-fold symmtery is formed by twisting one layer by \(30^\circ\) with respect to the center of \(\mathbf{c} = \frac{2}{3}\mathbf{a}_1 + \frac{2}{3}\mathbf{a}_2\), where \(\mathbf{a}_1\) and \(\mathbf{a}_2\) are the lattice vectors of the primitive cell of fixed layer. The script can be found at examples/advanced/quasi_crystal.py. We begin with importing the packages and defining the geometric parameters:

import math

import numpy as np
from numpy.linalg import norm

import tbplas as tb


angle = 30 / 180 * math.pi
center = (2./3, 2./3, 0)
radius = 3.0
shift = 0.3349
dim = (33, 33, 1)

Here angle is the twisting angle and center is the fractional coordinate of twisting center. The radius of the quasicrystal is controlled by radius, while shift specifies the interlayer distance. We need a large cell to hold the quasicrystal, whose dimension is given in dim. After introducing the parameters, we firstly get the Cartesian coordinate of the twisting center by:

# Get the Cartesian coordinate of twisting center
prim_cell = tb.make_graphene_diamond()
center = np.array([dim[0]//2, dim[1]//2, 0]) + center
center = np.matmul(center, prim_cell.lat_vec)

Since we have extended the primitive cell by \(33\times33\times1\) times, and we want the quasicrystal to be located in the center of the cell, we need to convert the coordinate of twisting center in line 2-3. Then we build the fixed and twisted layers by:

# Build fixed and twisted layers
layer_fixed = tb.extend_prim_cell(prim_cell, dim=dim)
layer_twisted = tb.extend_prim_cell(prim_cell, dim=dim)

and shift and twist the twisted layer with respect to the center:

# Twist and shift twisted layer
tb.spiral_prim_cell(layer_twisted, angle=angle, center=center, shift=shift)

The twisting operation is done by the spiral_prim_cell() function, where the Cartesian coordinate of the center is given in the center argument. Then we remove unnecessary orbitals to produce a round quasicrystal with finite radius. This is done by calling the cutoff_pc function:

# Remove unnecessary orbitals
cutoff_pc(layer_fixed, center=center, radius=radius)
cutoff_pc(layer_twisted, center=center, radius=radius)

which is defined as:

def cutoff_pc(prim_cell: tb.PrimitiveCell, center: np.ndarray,
              radius: float = 3.0) -> None:
    """
    Cutoff primitive cell up to given radius with respect to given center.

    :param prim_cell: supercell to cut
    :param center: Cartesian coordinate of center in nm
    :param radius: cutoff radius in nm
    :return: None. The incoming supercell is modified.
    """
    idx_remove = []
    orb_pos = prim_cell.orb_pos_nm
    for i, pos in enumerate(orb_pos):
        if norm(pos[:2] - center[:2]) > radius:
            idx_remove.append(i)
    prim_cell.remove_orbitals(idx_remove)
    prim_cell.trim()

where we loop over orbital positions to collect the indices of unnecessary orbitals, then call the remove_orbitals and trim functions. Before merging the fixed and twisted layers, we need to reset the lattice vectors and origin of twisted layer to that of fixed layer by calling the reset_lattice function:

# Reset the lattice of twisted layer
layer_twisted.reset_lattice(layer_fixed.lat_vec, layer_fixed.origin,
                            unit=tb.NM, fix_orb=True)

After that, we can merge them safely:

# Merge layers
final_cell = tb.merge_prim_cell(layer_fixed, layer_twisted)

Finally, we extend the hoppings and visualize the quasicrystal:

# Extend and visualize the model
extend_hop(final_cell)
final_cell.plot(with_cells=False, with_orbitals=False, hop_as_arrows=False, hop_eng_cutoff=0.3)

where the extend_hop function is defined in Build hetero-structure. The output is shown in following figure:

Plot of the quasicrystal formed from the incommensurate \(30^\circ\) twisted bilayer graphene with a radius of 3 nm.