Effective Dynamics of 2D Bloch Electrons in External Fields Derived From Symmetry

TL;DR

This paper develops a symmetry-based framework for deriving effective dynamics of 2D Bloch electrons, including Hamiltonians that account for spin-orbit coupling, strain, and external fields, with applications to MoS2 and graphene.

Contribution

It introduces a systematic symmetry-based method to derive multiband Hamiltonians for 2D materials, accounting for coordinate system choices and external influences.

Findings

Symmetry-adapted basis functions block-diagonalize the TB Hamiltonian.

Coordinate system choices affect IRs at the Brillouin zone boundary.

Constructed a multiband Hamiltonian for MoS2 near the K point.

Abstract

We develop a comprehensive theory for the effective dynamics of Bloch electrons based on symmetry. We begin with a scheme to systematically derive the irreducible representations (IRs) characterizing the Bloch functions. Starting from a tight-binding (TB) approach, we decompose the TB basis functions into localized symmetry-adapted atomic orbitals and crystal-periodic symmetry-adapted plane waves. Each of these subproblems is independent of the details of a particular crystal structure and it is largely independent of the other subproblem, hence permitting for each subproblem an independent universal solution. Taking monolayer MoS$_2$ and few-layer graphene as examples, we tabulate the symmetrized $p$ and $d$ orbitals as well as the symmetrized plane wave spinors for these systems. The symmetry-adapted basis functions block-diagonalize the TB Hamiltonian such that each block yields eigenstates transforming according to one of the IRs of the group of the wave vector $G_k$. For many crystal structures, it is possible to define multiple distinct coordinate systems such that for wave vectors $k$ at the border of the Brillouin zone the IRs characterizing the Bloch states depend on the coordinate system, i.e., these IRs of $G_k$ are not uniquely determined by the symmetry of a crystal structure. The different coordinate systems are related by a coordinate shift that results in a rearrangement of the IRs of $G_k$ characterizing the Bloch states. We illustrate this rearrangement with three coordinate systems for MoS$_2$ and tri-layer graphene. Using monolayer MoS$_2$ as an example, we combine the symmetry analysis of its bulk Bloch states with the theory of invariants to construct a generic multiband Hamiltonian for electrons near the $K$ point of the Brillouin zone. The Hamiltonian includes the effect of spin-orbit coupling, strain and external electric and magnetic fields.

Equations

$$g\exp (i\mathrm{k}\cdot \mathrm{r})=\exp (i\mathrm{k}\cdot \mathrm{g}\cdot \mathrm{r})\equiv \exp (i\mathrm{k}\cdot {\mathrm{r}}^{\prime })$$

$${\mathrm{r}}^{\prime }=\mathrm{g}\cdot \mathrm{r}.$$

$$g\exp (i\mathrm{k}\cdot \mathrm{r})$$

$$=\exp (i{\mathrm{k}}^{\prime }\cdot \mathrm{r})$$

$${\mathrm{k}}^{\prime }={\mathrm{g}}^{-1}\cdot \mathrm{k}.$$

$${\mathrm{k}}^{\prime }={\mathrm{g}}^{-1}\cdot \mathrm{k}=\mathrm{k}+{\mathrm{b}}_g,$$

$$g\exp (i\mathrm{k}\cdot \mathrm{a})=\exp (i{\mathrm{k}}^{\prime }\cdot \mathrm{a})=\exp (i\mathrm{k}\cdot \mathrm{a})$$

$${\Phi }_{\nu \mathrm{k}}^{W\mu }(\mathrm{r})=\frac{e^{i\mathrm{k}\cdot \mathrm{r}}}{\sqrt{N}}\sum _j{e}^{-i\mathrm{k}\cdot (\mathrm{r}-{\mathrm{R}}_j^{W\mu })}{\varphi }_{\nu }^W(\mathrm{r}-{\mathrm{R}}_j^{W\mu }),$$

$$H(\mathrm{k}{)}_{\nu {\nu }^{\prime }}^{W\mu W^{\prime }{\mu }^{\prime }}$$

$$={ϵ}_{\nu }^W{\delta }_{\nu {\nu }^{\prime }}{\delta }_{W{W}^{\prime }}{\delta }_{\mu {\mu }^{\prime }}+\text{sideset}{}^{\prime }{\sum }_{j{j}^{\prime }}{t}_{\nu {\nu }^{\prime }j{j}^{\prime }}^{W{W}^{\prime }\mu {\mu }^{\prime }},$$

$${ϵ}_{\nu }^W$$

$$=\int {\varphi }_{\nu }^{W\ast }(\mathrm{r})H{\varphi }_{\nu }^W(\mathrm{r})d^{3}r$$

$$t_{\nu {\nu }^{\prime }j{j}^{\prime }}^{W{W}^{\prime }\mu {\mu }^{\prime }}$$

$$\times \int {\varphi }_{\nu }^{W\ast }(\mathrm{r}-{\mathrm{R}}_j^{W\mu })H{\varphi }_{{\nu }^{\prime }}^{W^{\prime }}(\mathrm{r}-{\mathrm{R}}_{j^{\prime }}^{W^{\prime }{\mu }^{\prime }})d^{3}r$$

$${\tilde{\varphi }}_{\nu \mathrm{k}}^W(\mathrm{r}-{\mathrm{R}}_j^{W\mu })\equiv e^{-i\mathrm{k}\cdot (\mathrm{r}-{\mathrm{R}}_j^{W\mu })}{\varphi }_{\nu }^W(\mathrm{r}-{\mathrm{R}}_j^{W\mu }),$$

$${\Phi }_{\nu \mathrm{k}}^{W\mu }(\mathrm{r})=\frac{e^{i\mathrm{k}\cdot \mathrm{r}}}{\sqrt{N}}\sum _j{\tilde{\varphi }}_{\nu \mathrm{k}}^W(\mathrm{r}-{\mathrm{R}}_j^{W\mu }).$$

$${\tilde{\varphi }}_{\nu \mathrm{k}}^W(\mathrm{r}-{\mathrm{R}}_j^{W\mu })\approx {\varphi }_{\nu }^W(\mathrm{r}-{\mathrm{R}}_j^{W\mu }).$$

$${\Phi }_{\nu \mathrm{k}}^{W\mu }(\mathrm{r})$$

$$=e^{i\mathrm{k}\cdot \mathrm{r}}{\mathcal{A}}_{\nu }^{W\mu }(\mathrm{r})$$

$${\mathcal{A}}_{\nu }^{W\mu }(\mathrm{r})=\frac{1}{\sqrt{N}}\sum _j{\varphi }_{\nu }^W(\mathrm{r}-{\mathrm{R}}_j^{W\mu })$$

$${\mathcal{A}}_{\nu }^{W\mu }(\mathrm{r}={\mathrm{R}}_j^{W\mu }+\delta \mathrm{r})=\frac{1}{\sqrt{N}}\sum _j{\varphi }_{\nu }^W(\delta \mathrm{r})\propto {\varphi }_{\nu }^W(\delta \mathrm{r}),$$

$$\exp [i\mathrm{k}\cdot ({\mathrm{R}}_j^{W\mu }+\delta \mathrm{r})]\approx \exp (i\mathrm{k}\cdot {\mathrm{R}}_j^{W\mu })\equiv q_{\mathrm{k}}({\mathrm{R}}_j^{W\mu }),$$

$${\Phi }_{\nu \mathrm{k}}^{W\mu }(\mathrm{r}={\mathrm{R}}_j^{W\mu }+\delta \mathrm{r})\approx q_{\mathrm{k}}({\mathrm{R}}_j^{W\mu }){\varphi }_{\nu }^W(\delta \mathrm{r}).$$

$${\Psi }_{n\mathrm{k}}(\mathrm{r})=\sum _{W\mu \nu }{\psi }_{\nu n\mathrm{k}}^{W\mu }{\Phi }_{\nu \mathrm{k}}^{W\mu }(\mathrm{r})$$

$${\Psi }_{n\mathrm{k}}(\mathrm{r}={\mathrm{R}}_j^{W\mu }+\delta \mathrm{r})\approx \sum _{W\mu }{q}_{\mathrm{k}}({\mathrm{R}}_j^{W\mu })\sum _{\nu }{\psi }_{\nu n\mathrm{k}}^{W\mu }{\varphi }_{\nu }^W(\delta \mathrm{r}).$$

$${\Gamma }_{\mathrm{k}W}^{\Phi }\equiv {\Gamma }_{\mathrm{k}W}^q\times {\Gamma }_{\mathrm{k}W}^{\varphi }.$$

$${\Gamma }_{\mathrm{k}W}^{\Phi }=\sum _I{a}_I{\Gamma }_I,$$

$${\Phi }_{\mathrm{k}}^{WI\alpha \beta }(\mathrm{r})=\sum _{\mu \nu }{\mathcal{C}}_{I\alpha \beta }^{W\mu \nu \mathrm{k}}{\Phi }_{\nu \mathrm{k}}^{W\mu }(\mathrm{r}),$$

$$H({\bm{\mathrm{k}}})_{\alpha\alpha^{\prime}\beta\beta^{\prime}}^{WW^{\prime}II^{\prime}}=\delta_{II^{\prime}}\Bigl{[}\epsilon_{I\alpha}^{W}\,\delta_{WW^{\prime}}\,\delta_{\alpha\alpha^{\prime}}\,\delta_{\beta\beta^{\prime}}+\sideset{}{{}^{\prime}}{\sum}_{jj^{\prime}}t_{\beta\beta^{\prime}\,jj^{\prime}}^{WW^{\prime}I\,\alpha\alpha^{\prime}}\Bigr{]}.$$

$${ϵ}_{I\alpha }^W\equiv \int {\Phi }_{\mathrm{k}}^{WI\alpha \beta \ast }(\mathrm{r})H{\Phi }_{\mathrm{k}}^{WI\alpha \beta }(\mathrm{r})d^{3}r$$

$$t_{\beta {\beta }^{\prime }j{j}^{\prime }}^{W{W}^{\prime }I\alpha {\alpha }^{\prime }}\equiv \int {\Phi }_{\mathrm{k}}^{WI\alpha \beta \ast }(\mathrm{r})H{\Phi }_{\mathrm{k}}^{W^{\prime }I{\alpha }^{\prime }{\beta }^{\prime }}({\mathrm{r}}^{\prime })d^{3}r$$

$${\Psi }_{n\mathrm{k}}^{I\beta }(\mathrm{r})=\sum _{W,\alpha }{C}_{n\mathrm{k}}^{WI\alpha }{\Phi }_{\mathrm{k}}^{WI\alpha \beta }(\mathrm{r})$$

$$g{\varphi }_{\nu }^W(\delta \mathrm{r})={\varphi }_{\nu }^W(g\delta \mathrm{r})=\sum _{{\nu }^{\prime }}{\mathcal{D}}_{\mathrm{k}W}^{\varphi }(g{)}_{\nu {\nu }^{\prime }}{\varphi }_{{\nu }^{\prime }}^W(\delta \mathrm{r}),$$

$${\varphi }_{I\beta }^W(\delta \mathrm{r})$$

Full text

For many crystal structures, it is possible to define multiple distinct coordinate systems such that for wave vectors ${\bm{\mathrm{k}}}$ at the border of the Brillouin zone the IRs characterizing the Bloch states depend on the coordinate system, i.e., these IRs of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ are not uniquely determined by the symmetry of a crystal structure. The different coordinate systems are related by a coordinate shift that results in a rearrangement of the IRs of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ characterizing the Bloch states. We illustrate this rearrangement with three coordinate systems for MoS2 and trilayer graphene.

The freedom to choose different distinct coordinate systems can simplify the symmetry analysis of the Bloch states. Given the IRs of the Bloch states in one coordinate system, a rearrangement lemma yields immediately the IRs of the Bloch states in the other coordinate systems. The rearrangement of the IRs in different coordinate systems does not affect observable physics such as selection rules or the effective Hamiltonians for the dynamics of Bloch states in external fields.

Using monolayer MoS2 as an example, we combine the symmetry analysis of its bulk Bloch states with the theory of invariants to construct a generic multiband Hamiltonian for electrons near the ${\bm{\mathrm{K}}}$ point of the Brillouin zone. The Hamiltonian includes the effect of spin-orbit coupling, strain and external electric and magnetic fields. Invariance of the Hamiltonian under time reversal yields additional constraints for the allowed terms in the Hamiltonian and it determines the phase (real or imaginary) of the prefactors.

I Introduction

Near a band extremum, the electron dynamics in a crystalline solid resembles the dynamics of free electrons in the absence of the periodic crystal potential. In the multiband envelope-function approximation (EFA) the electrons are characterized by an $N\times N$ Hamiltonian $\mathcal{H}$ for $N$-component spinors conceptually similar to relativistic electrons described by the Dirac equation [1, 2, 3, 4, 5]. The simplest approach within the EFA is the effective-mass approximation (EMA) that represents the electron dynamics by a Schrödinger equation with effective mass $m^{\ast}$ reflecting the curvature of the band dispersion $E({\bm{\mathrm{k}}})$. External electric and magnetic fields ${\bm{\mathcal{E}}}$ and ${\bm{\mathcal{B}}}$ break the lattice periodicity of the crystal structure. It is an important advantage of EFA and EMA that they allow one to incorporate the field ${\bm{\mathcal{E}}}$ by adding the corresponding scalar potential $\Phi$ to the diagonal of the Hamiltonian, and the operator of crystal momentum $\hbar{\bm{\mathrm{k}}}=-i\hbar\nabla$ is replaced by $-i\hbar\nabla+e{\bm{\mathrm{A}}}$, where ${\bm{\mathrm{A}}}$ is the vector potential for the magnetic field ${\bm{\mathcal{B}}}$. Other perturbations such as spin-orbit coupling, strain and electron-phonon coupling can likewise be included in the Hamiltonian [3]. This is a major reason why EFA and EMA are very popular for theoretical studies of both bulk semiconductors (e.g., Refs. [6, 7, 3, 8, 9, 10, 11]) and semiconductor quantum structures (e.g., Refs. [4, 12, 5, 13, 14, 15, 16]).

The form of the Hamiltonian $\mathcal{H}$ depends on the symmetry of the crystal structure and more specifically on the symmetry of the bulk electronic states that are included in $\mathcal{H}$ [17, 18, 3]. The relevant symmetry group for states with wave vector ${\bm{\mathrm{k}}}$ is the point group $\mathcal{G}_{{\bm{\mathrm{k}}}}$ which includes those symmetry elements of the crystallographic point group $\mathcal{G}_{0}$ (crystal class) which either leave ${\bm{\mathrm{k}}}$ unchanged or map ${\bm{\mathrm{k}}}$ onto an equivalent wave vector. The symmetry of individual states at ${\bm{\mathrm{k}}}$ is characterized by the respective irreducible representations (IRs) of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ according to which these states transform. The general form of the Hamiltonian $\mathcal{H}$ including its dependence on, e.g., spin-orbit coupling, strain and external fields can then be derived from its invariance under $\mathcal{G}_{{\bm{\mathrm{k}}}}$ [17, 3]. Here we develop a general theory to determine the IRs of Bloch functions for a given wave vector ${\bm{\mathrm{k}}}$, focusing for clarity on symmorphic space groups. Using a tight-binding (TB) approach along with the fact that the atomic orbitals are localized in the vicinity of the atomic sites we demonstrate that the TB basis functions can be factorized into localized symmetry-adapted atomic orbitals and crystal-periodic symmetry-adapted plane waves. Each of these two subproblems permits a universal classification, independent of the details of a particular crystal structure and also largely independent of the other subproblem. The symmetrized atomic orbitals depend only on the angular momentum of the atomic orbitals and the point group $\mathcal{G}_{{\bm{\mathrm{k}}}}$ of the wave vector ${\bm{\mathrm{k}}}$; but these orbitals are independent of the specific type of atom and the details of the crystal structure. The symmetrized plane waves form discrete Bloch functions that depend on the wave vector ${\bm{\mathrm{k}}}$ and the Wyckoff positions of the atoms in a crystal structure; but they are independent of the type of atoms occupying these positions. The symmetry-adapted basis functions block-diagonalize the TB Hamiltonian such that each block yields eigenstates transforming according to one of the IRs of the group of the wave vector $\mathcal{G}_{{\bm{\mathrm{k}}}}$.

Given the symmetry group $G$ of a quantum system, the IRs of $G$ are generally assumed to provide a distinct label for the eigenstates of the Hamiltonian, as noted by Wigner: “The representation of the group of the Schrödinger equation which belongs to a particular eigenvalue is uniquely determined up to a similarity transformation.” (Ref. [19], p. 110, highlighting adopted from Ref. [19]). This uniqueness of the IRs is immediately relevant for many physical properties of a physical system that depend on the symmetry of its electronic states. For example, the Wigner-Eckart theorem allows one to express the selection rules for optical transitions in terms of the IRs of the initial and final states between which a transition occurs [19]. Similarly, the EFA Hamiltonians $\mathcal{H}$ depend on the IRs of the bands described by $\mathcal{H}$, as noted above. We demonstrate that the IRs characterizing the Bloch eigenstates in certain crystals including transition metal dichalcogenides (TMDCs) are not unique, but they depend on the coordinate system used to describe the space group symmetries of these materials [20, 21, 22, 23]. We show that distinct valid coordinate systems are related by a coordinate shift that defines a rearrangement representation. The IRs of the electronic states in the different coordinate systems are then related via a rearrangement lemma that facilitates the symmetry analysis of Bloch states. Also, we show how important physics including optical selection rules and EFA Hamiltonians $\mathcal{H}$, despite the rearrangement of band IRs, does not depend on the coordinate system being used.

Our general theory applies to any crystalline material. For a detailed example, we focus on a monolayer of the TMDC MoS2. TMDCs are of the general form $TX_{2}$, where $T$ is a transition-metal such as Mo or W and $X$ is a chalcogen which can be S, Se, or Te. Three-dimensional (3D) bulk $TX_{2}$ consists of covalently bonded 2D monolayers coupled vertically by weak van der Waals forces [24], making it possible to obtain monolayers via, e.g., mechanical exfoliation [25]. Electronic band structure calculations have shown that bulk 2H-MoS2 is a semiconductor [24]. More recently, optical spectroscopy [25] and theoretical studies [26, 27, 28] found that decreasing the number of layers changes the fundamental gap from indirect to direct in the limit of a single monolayer. The spin-dependent dispersion of monolayer TMDCs has been studied using TB [29, 30] and ${\bm{\mathrm{k}}}\cdot{\bm{\mathrm{p}}}$ methods [31, 16]. See Refs. [32, 33, 34, 35] for general reviews of 2D TMDCs. In this paper, we combine our symmetry analysis for the bulk Bloch states in monolayer MoS2 with the theory of invariants [3] to derive a generic multiband EFA Hamiltonian for electrons near the ${\bm{\mathrm{K}}}$ point of the Brillouin zone (BZ). The Hamiltonian includes the effect of strain, external electric and magnetic fields, spin and valley degrees of freedom. For comparison, we also perform a symmetry analysis for few-layer graphene [14] which confirms earlier work [36, 37]. We note that our work expands on the theory of IRs for point and space groups [38, 3, 39]. It is conceptually rather different from recent work on band representations [40, 41, 42].

In Sec. II, we develop the general theory of the symmetry of TB Bloch functions. The decomposition of TB wave functions is discussed in Sec. II.3 followed by detailed discussions of the symmetrized atomic orbitals (Sec. II.5) and symmetrized plane waves (Sec. II.6). The rearrangement of the IRs of Bloch states under a change of the coordinate system is discussed in Sec. II.7. We use the general formalism of Sec. II to derive the symmetry of bulk Bloch states in monolayer TMDCs (Sec. III) such as MoS2 and to few-layer graphene (Sec. IV). We show in Sec. V how optical selection rules are not affected by the rearrangement of IRs under a change of coordinate system. In Sec. VI we derive the generic invariant expansion of the EFA Hamiltonian $\mathcal{H}$ for MoS2. Section VII contains our conclusions.

II Symmetry of Bloch Functions

Very generally, the eigenstates of a Hamiltonian transform according to an IR of the symmetry group of the Hamiltonian. In band theory it is thus an important goal to determine the IRs of the energy bands $E_{n}({\bm{\mathrm{k}}})$ and corresponding Bloch functions $\Psi_{n{\bm{\mathrm{k}}}}({\bm{\mathrm{r}}})$, where the symmetry group $\mathcal{G}_{{\bm{\mathrm{k}}}}$ at a given wave vector ${\bm{\mathrm{k}}}$ is called the point group of the wave vector ${\bm{\mathrm{k}}}$. In this section we discuss a general method for determining the transformation properties of Bloch functions with a certain wave vector ${\bm{\mathrm{k}}}$, which allows us to determine the corresponding IRs of $\mathcal{G}_{{\bm{\mathrm{k}}}}$. We also discuss a rearrangement lemma for the IRs characterizing the Bloch functions in a crystal. Applications to specific materials such as monolayer MoS2 will be discussed in subsequent sections.

II.1 The group of the wave vector

In the following, we will repeatedly need to evaluate the action of a point symmetry operation $g$ on a plane wave $\exp(i{\bm{\mathrm{k}}}\cdot{\bm{\mathrm{r}}})$. Here, $g$ can be represented via an orthogonal $3\times 3$ matrix ${\bm{\mathrm{g}}}$. (In the context of quasi-2D materials discussed below ${\bm{\mathrm{g}}}$ becomes a $2\times 2$ matrix.) Thus we have

[TABLE]

Note that when transforming the position vector ${\bm{\mathrm{r}}}$, the wave vector ${\bm{\mathrm{k}}}$ is a fixed parameter characterizing the plane wave $\exp(i{\bm{\mathrm{k}}}\cdot{\bm{\mathrm{r}}})$ that does not change under $g$. Nonetheless, since $g$ is an orthogonal transformation, we can also write Eq. (1a) as

[TABLE]

Thus we can evaluate $g$ either by transforming the position vector ${\bm{\mathrm{r}}}$ or by inversely transforming the wave vector ${\bm{\mathrm{k}}}$.

In the group theory of crystallographic space groups, the point-group symmetries $g$ of Bloch functions $\Psi_{\bm{\mathrm{k}}}({\bm{\mathrm{r}}})$ with wave vector ${\bm{\mathrm{k}}}$ form the point group $\mathcal{G}_{{\bm{\mathrm{k}}}}$ of the wave vector ${\bm{\mathrm{k}}}$ [38, 3, 39]. Given the point group $\mathcal{G}_{0}$ of a crystal structure, the group $\mathcal{G}_{{\bm{\mathrm{k}}}}$ is defined by the condition that it contains the symmetry elements of $\mathcal{G}_{0}$ that map ${\bm{\mathrm{k}}}$ onto a vector $\bm{\mathrm{k}}^{\prime}$ such that

[TABLE]

where ${\bm{\mathrm{b}}}_{g}$ is a reciprocal lattice vector with the possibility ${\bm{\mathrm{b}}}_{g}=0$. Indeed, since $g$ represents point group operations, we can have ${\bm{\mathrm{b}}}_{g}\neq 0$ only if ${\bm{\mathrm{k}}}$ is from the border of the BZ. For positions ${\bm{\mathrm{r}}}={\bm{\mathrm{a}}}$ that are lattice vectors we have

[TABLE]

by definition of $\mathcal{G}_{{\bm{\mathrm{k}}}}$.

II.2 Tight-binding Hamiltonian

We denote the TB basis functions (that are Bloch functions) as

[TABLE]

where $\phi_{\nu}^{W}({\bm{\mathrm{r}}}-{\bm{\mathrm{R}}}_{j}^{W\mu})$ are the atomic orbitals of type $\nu$ centered about the positions ${\bm{\mathrm{R}}}_{j}^{W\mu}$ of the atoms. The label $W$ denotes the Wyckoff letter of the atomic positions of the crystal structure [43]. The label $\mu$ has values $\mu=1,\dots,m$, where $m$ is the multiplicity of $W$. The index $j$ labels the unit cells of the crystal structure; it runs through the positions in a Bravais lattice. The matrix elements of the TB Hamiltonian can then be written as

[TABLE]

where

[TABLE]

denotes the on-site energies for the atomic orbitals (that do not depend on the indices $\mu$ and $j$) and

[TABLE]

are the hopping integrals. The prime on the summation sign in Eq. (6b) indicates that the sum excludes the on-site term $\epsilon_{\nu}^{W}$. The TB approximation implies that this sum is restricted to $n$th-nearest neighbors with a small value of $n$. The hopping integrals can be written in terms of the Slater-Koster parameters [44] for hopping integrals between atomic orbitals at positions ${\bm{\mathrm{R}}}_{j}^{W\mu}$ and ${\bm{\mathrm{R}}}_{j^{\prime}}^{W^{\prime}\mu^{\prime}}$.

II.3 Decomposition of TB wave functions

Generally, the basis functions (5) for a given wave vector ${\bm{\mathrm{k}}}$ transform according to a representation $\Gamma^{\Phi}_{{\bm{\mathrm{k}}}W}$ of the group of the wave vector $\mathcal{G}_{{\bm{\mathrm{k}}}}$ that need not be irreducible. (Here the generic superscript $\Phi$ accounts for the fact that multiple atomic orbitals with different indices $\nu$ may transform jointly according to the same representation. By definition of Wyckoff letters $W$, orbitals at different positions $\mu$ of a given Wyckoff letter $W$ transform jointly according to one representation.) As a consequence of the Wigner-Eckart theorem [19, 45] and the fact that $H$ transforms according to the identity representation $\Gamma_{1}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$, the hopping integrals (8) vanish when the product of the representations $\Gamma^{\Phi}_{{\bm{\mathrm{k}}}W}$ and $\Gamma^{\Phi^{\prime}}_{{\bm{\mathrm{k}}}W^{\prime}}$ does not contain the identity representation. In the following, we present a general scheme for transforming the set of basis functions (5) into a symmetry-adapted set of basis functions, where each function transforms irreducibly under $\mathcal{G}_{{\bm{\mathrm{k}}}}$, so that two such basis functions only couple when they both transform according to the same IR of $\mathcal{G}_{{\bm{\mathrm{k}}}}$. This scheme is based on a decomposition of the basis functions into symmetry-adapted plane waves and symmetry-adapted atomic orbitals.

We denote

[TABLE]

so that Eq. (5) becomes

[TABLE]

Assuming for conceptual simplicity that the atomic orbitals $\phi_{\nu}^{W}$ are localized over a region much smaller than the nearest-neighbor distance 111This assumption facilitates a discussion of the symmetry of the basis functions (5). Functions (5) with a finite overlap between nearest neighbors that are needed in a TB calculation must have the same symmetry as the simplified functions discussed here. The symmetry of these functions cannot change discontinuously when the overlap between the atomic orbitals is switched off., the functions $\tilde{\phi}_{\nu{\bm{\mathrm{k}}}}^{W}({\bm{\mathrm{r}}}-{\bm{\mathrm{R}}}_{j}^{W\mu})$ are only nonzero for ${\bm{\mathrm{r}}}$ close to ${\bm{\mathrm{R}}}_{j}^{W\mu}$. Hence, in the vicinity of any atomic position ${\bm{\mathrm{R}}}_{j}^{W\mu}$, i.e., for ${\bm{\mathrm{r}}}\equiv{\bm{\mathrm{R}}}_{j}^{W\mu}+\delta{\bm{\mathrm{r}}}$ with small $\delta{\bm{\mathrm{r}}}$, we have

[TABLE]

Therefore, the TB basis function $\Phi_{\nu{\bm{\mathrm{k}}}}^{W\mu}({\bm{\mathrm{r}}})$ can be approximated as

[TABLE]

with atomic functions

[TABLE]

independent of the wave vector ${\bm{\mathrm{k}}}$. For strongly localized atomic orbitals and positions ${\bm{\mathrm{r}}}={\bm{\mathrm{R}}}_{j}^{W\mu}+\delta{\bm{\mathrm{r}}}$ we have

[TABLE]

that is, near an atomic site ${\bm{\mathrm{R}}}_{j}^{W\mu}$, the atomic functions $\mathcal{A}_{\nu}^{W\mu}({\bm{\mathrm{r}}})$ are simply proportional to $\phi_{\nu}^{W}(\delta{\bm{\mathrm{r}}})$, independent of the index $\mu$. Therefore, the atomic functions $\mathcal{A}_{\nu}^{W\mu}({\bm{\mathrm{r}}}={\bm{\mathrm{R}}}_{j}^{W\mu}+\delta{\bm{\mathrm{r}}})$ have the same symmetry properties as the atomic orbitals $\phi_{\nu}^{W}(\delta{\bm{\mathrm{r}}})$.

The plane wave $e^{i{\bm{\mathrm{k}}}\cdot{\bm{\mathrm{r}}}}$ for positions ${\bm{\mathrm{r}}}={\bm{\mathrm{R}}}_{j}^{W\mu}+\delta{\bm{\mathrm{r}}}$ near an atomic site ${\bm{\mathrm{R}}}_{j}^{W\mu}$ is approximately given by

[TABLE]

where $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ denotes the plane wave with wave vector ${\bm{\mathrm{k}}}$ associated with the Wyckoff position ${\bm{\mathrm{R}}}_{j}^{W\mu}$ for fixed $W$ and $\mu$, but $j$ runs over all positions in a Bravais lattice. These discrete quantities $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ will be discussed in more detail in Sec. II.6. The TB basis function thus can be factorized (ignoring normalization)

[TABLE]

This expression will be analyzed further in the following sections.

As a side remark, we note that the eigenfunctions of the TB Hamiltonian (6) for the band $n$ and wave vector ${\bm{\mathrm{k}}}$ expressed in terms of the basis functions (5) take the form

[TABLE]

with expansion coefficients $\psi_{\nu n{\bm{\mathrm{k}}}}^{W\mu}$. These eigenfunctions permit a factorization similar to Eq. (16)

[TABLE]

However, the discussion of the TB wave functions $\Psi_{n{\bm{\mathrm{k}}}}({\bm{\mathrm{r}}})$ is greatly simplified if instead of the basis functions (5) we use symmetry-adapted basis functions to be discussed in the following.

II.4 Symmetry-adapted basis functions

The main advantage of the approximate expression (16) lies in the fact that the function (16) has the same symmetry properties as the TB basis function (5). Yet the factorization in Eq. (16) allows one to discuss the symmetry of the plane waves $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ (characterized by a representation $\Gamma^{q}_{{\bm{\mathrm{k}}}W}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$, see Sec. II.6) separate from the symmetry of the atomic orbitals $\phi_{\nu}^{W}(\delta{\bm{\mathrm{r}}})$ (characterized by a representation $\Gamma^{\phi}_{{\bm{\mathrm{k}}}W}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$, see Sec. II.5). Often, these representations are irreducible, though generally they can be reducible. The TB basis functions (5) at one Wyckoff letter $W$ then transform according to the product representation

[TABLE]

The representation $\Gamma^{\Phi}_{{\bm{\mathrm{k}}}W}$ may be reducible (even if $\Gamma^{q}_{{\bm{\mathrm{k}}}W}$ and $\Gamma^{\phi}_{{\bm{\mathrm{k}}}W}$ are irreducible), giving the decomposition

[TABLE]

where $a_{I}\geq 0$ are the multiplicities with which the IRs $\Gamma_{I}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ appear in $\Gamma^{\Phi}_{{\bm{\mathrm{k}}}W}$. Generally, the IRs $\Gamma_{I}$ contained in $\Gamma^{\Phi}_{{\bm{\mathrm{k}}}W}$ give symmetry-adapted basis functions

[TABLE]

where $\alpha=1,2,\ldots,a_{I}$ and $\beta$ labels the different functions transforming jointly according to $\Gamma_{I}$. The coefficients $\mathcal{C}_{I\alpha\beta}^{W\mu\nu{\bm{\mathrm{k}}}}$ represent the weights with which the functions $\Phi_{\nu{\bm{\mathrm{k}}}}^{W\mu}({\bm{\mathrm{r}}})$ contribute to the symmetry-adapted basis functions $\Phi_{\bm{\mathrm{k}}}^{WI\alpha\beta}({\bm{\mathrm{r}}})$. Given symmetry-adapted plane waves transforming according to the IR $\Gamma_{J}$ and atomic orbitals transforming according to the IR $\Gamma_{J^{\prime}}$ (see discussion below), the expansion coefficients $\mathcal{C}_{I\alpha\beta}^{W\mu\nu{\bm{\mathrm{k}}}}$ are given by the Clebsch-Gordan coefficients for coupling $\Gamma_{J}$ and $\Gamma_{J^{\prime}}$ to obtain $\Gamma_{I}$.

For the symmetry-adapted basis functions (21), the matrix elements (6) of the TB Hamiltonian become block-diagonal with respect to different IRs $\Gamma_{I}$

[TABLE]

Here

[TABLE]

denotes the on-site energies. These can always be made diagonal in the index $\alpha$ by a suitable definition of the $a_{I}$ sets of basis functions (21) transforming according to $\Gamma_{I}$. The hopping matrix elements become

[TABLE]

with ${\bm{\mathrm{r}}}={\bm{\mathrm{R}}}_{j}^{W\mu}+\delta{\bm{\mathrm{r}}}$ and ${\bm{\mathrm{r}}}^{\prime}={\bm{\mathrm{R}}}_{j^{\prime}}^{W^{\prime}\mu^{\prime}}+\delta{\bm{\mathrm{r}}}$. We remark that in actual TB models it may happen that a block (22) can be further decomposed into subblocks if symmetry-allowed couplings between distant neighbors are ignored within the TB approximation.

Each block (22) of the TB Hamiltonian yields eigenfunctions

[TABLE]

with expansion coefficients $C_{n{\bm{\mathrm{k}}}}^{WI\alpha\beta}$ that transform irreducibly according to the $\beta$th component of the IR $\Gamma_{I}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$. To proceed, we discuss first the symmetry of the atomic orbitals $\phi_{\nu}^{W}(\delta{\bm{\mathrm{r}}})$, followed by a discussion of the symmetry of the plane waves $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$.

II.5 Transformation of atomic orbitals $\phi_{\nu}^{W}(\delta\mathbf{r})$

We can study the symmetry of the atomic function $\mathcal{A}_{\nu}^{W\mu}({\bm{\mathrm{r}}})$ in the vicinity of the atomic positions ${\bm{\mathrm{R}}}_{j}^{W\mu}$ by only looking at the orbitals $\phi_{\nu}^{W}(\delta{\bm{\mathrm{r}}})$. Obviously, this problem is independent of the index $\mu$, the band index $n$, and the wave vector ${\bm{\mathrm{k}}}$ though, of course, we are generally interested in the transformational behavior of these orbitals with respect to the group $\mathcal{G}_{{\bm{\mathrm{k}}}}$ of the wave vector ${\bm{\mathrm{k}}}$. For symmetry operations $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$ we have

[TABLE]

where $\Gamma^{\phi}_{{\bm{\mathrm{k}}}W}=\{\mathcal{D}^{\phi}_{{\bm{\mathrm{k}}}W}(g):g\in\mathcal{G}_{{\bm{\mathrm{k}}}}\}$ is the representation describing the transformation of the atomic orbitals $\phi_{\nu}^{W}$. As usual, we assume that the atomic orbitals are characterized by some orbital angular momentum $l=0,1,2,\ldots$, so that the matrices $\mathcal{D}^{\phi}_{{\bm{\mathrm{k}}}W}(g)$ acquire a block structure corresponding to different values of $l$, indicating that there exists no mixing between atomic orbitals with different angular momenta under symmetry transformations $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$.

In general, the representation $\Gamma^{\phi}_{{\bm{\mathrm{k}}}W}$ is reducible 222The $2l+1$ atomic orbitals for a given magnitude $l$ of orbital angular momentum transform according to an IR of the full rotation group. The compatibility tables for decomposing the IRs of the full rotation group into IRs of the crystallographic point groups are listed in Ref. [56]. Here we require a more complete solution of this problem where we also construct the linear combinations of angular-momentum eigenstates that transform irreducibly according to the different IRs of $\mathcal{G}_{{\bm{\mathrm{k}}}}$. This problem does not depend on, e.g., the principal quantum numbers of these orbitals. We note that for sufficiently large $l$ (for any $\mathcal{G}_{{\bm{\mathrm{k}}}}$ if $l\geq 6$) these $2l+1$ atomic orbitals contribute to all (even or odd if $\mathcal{G}_{{\bm{\mathrm{k}}}}$ includes inversion) IRs of $\mathcal{G}_{{\bm{\mathrm{k}}}}$, which limits possibilities to decompose a TB Hamiltonian into blocks describing parts of the full band structure.. The projection operators (122) yield symmetry-adapted atomic orbitals $\phi^{W}_{I\beta}(\delta{\bm{\mathrm{r}}})$ transforming like component $\beta$ of the IR $\Gamma_{I}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$,

[TABLE]

We note that this analysis applies to the spinless case when the angular part of the atomic orbitals is given by the usual spherical harmonics $Y_{l}^{m}$. It can likewise be used in the spin-dependent case when the angular part of the atomic orbitals is given by spin-angular functions and the projection operators $\Pi_{I\beta}$ project on the double-group representations of $\mathcal{G}_{{\bm{\mathrm{k}}}}$.

II.6 Transformation of plane waves $q_{\mathbf{k}}(\mathbf{R}_{j}^{W\mu})$

While the IRs of TB eigenfunctions (17) depend on the band index $n$, the symmetry of the plane waves $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ can be discussed independent of the index $n$. Very generally, for a given Wyckoff letter $W$, the positions ${\bm{\mathrm{R}}}_{j}^{W\mu}$ transform among themselves under the operations of the space group. Hence, using the matrix ${\bm{\mathrm{g}}}$ defined in Eq. (1), we have

[TABLE]

i.e., a symmetry transformation $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$ generally changes both $\mu$ and $j$. We rewrite this as

[TABLE]

where the term in square brackets is a lattice vector. Applying the operation $g$ to the plane wave $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ yields

[TABLE]

Hence, the symmetry operation $g$ generally maps the plane wave $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ onto $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu^{\prime}})$ multiplied by a phase factor.

To analyze the mappings (30) further, we interpret the discrete plane waves $\{q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu}):\mu=1,\ldots,m\}$ as basis vectors in an $m$-dimensional vector space. Relative to this basis, we can express arbitrary points as $m$-component spinors. It facilitates this analysis to introduce $m$-component base spinors $\{\mathcal{q}_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu}):\mu=1,\ldots,m\}$ with components $\nu$ equal to $\delta_{\mu\nu}$. Using these base spinors, a plane wave at positions ${\bm{\mathcal{R}}}_{j}^{W}\equiv\{{\bm{\mathrm{R}}}_{j}^{W\mu}:\mu=1,\ldots,m\}$ becomes

[TABLE]

Applying a symmetry transformation $g$ to ${Q}_{\bm{\mathrm{k}}}({\bm{\mathcal{R}}}_{j}^{W})$ gives

[TABLE]

where $\mathcal{D}^{q}_{{\bm{\mathrm{k}}}W}(g)$ is an $m\times m$ matrix corresponding to the group element $g$. Each row and each column of $\mathcal{D}^{q}_{{\bm{\mathrm{k}}}W}(g)$ has only one nonzero matrix element which, according to Eq. (30), is given by

[TABLE]

where the $g$ dependence on the right-hand side is given by Eq. (29).

We show in the following that $\Gamma^{q}_{{\bm{\mathrm{k}}}W}\equiv\{\mathcal{D}^{q}_{{\bm{\mathrm{k}}}W}(g):g\in\mathcal{G}_{{\bm{\mathrm{k}}}}\}$ defines an $m$-dimensional representation of $\mathcal{G}_{{\bm{\mathrm{k}}}}$, i.e., for any two group elements $g_{1},g_{2}\in\mathcal{G}_{{\bm{\mathrm{k}}}}$ we have

[TABLE]

Given a position ${\bm{\mathrm{R}}}_{j_{0}}^{W\mu_{0}}$ and using Eq. (29), we have

[TABLE]

so that the phase (33) for the group elements $g_{1}$ and $g_{2}$ becomes

[TABLE]

The product $g_{2}\,g_{1}$, i.e., the transformation $g_{1}$ followed by $g_{2}$, is also an element of $\mathcal{G}_{{\bm{\mathrm{k}}}}$, giving

[TABLE]

As ${\bm{\mathrm{R}}}_{j_{1}}^{W\mu_{1}}-{\bm{\mathrm{R}}}_{j_{0}}^{W\mu_{1}}$ is a lattice vector, it follows from Eq. (4)

[TABLE]

with

[TABLE]

This confirms Eq. (34).

We remark that for ${\bm{\mathrm{k}}}=0$ the representation $\Gamma^{q}_{{\bm{\mathrm{k}}}W}$ is known as permutation representation [48] or equivalence representation [39], where it characterizes the permutations of $m$ objects under the symmetry operations $g\in\mathcal{G}_{0}$.

Using the fact that $g$ is an orthogonal transformation, we can also write Eq. (33) as

[TABLE]

For ${\bm{\mathrm{k}}}$ inside the Brillouin zone, where ${\bm{\mathrm{g}}}^{-1}\cdot{\bm{\mathrm{k}}}={\bm{\mathrm{k}}}$ by definition of $\mathcal{G}_{{\bm{\mathrm{k}}}}$, the nonzero matrix elements of $\mathcal{D}^{q}_{{\bm{\mathrm{k}}}W}(g)$ therefore become

[TABLE]

We see that the representation is generally nontrivial for $m>1$. However, for $m=1$, we have ${\bm{\mathrm{R}}}_{j}^{W\mu}={\bm{\mathrm{R}}}_{j}^{W\mu^{\prime}}$, so that $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ transforms according to the identity representation $\Gamma_{1}$.

II.6.1 Wyckoff positions with multiplicity $m=1$

We consider first the special case of Wyckoff positions with multiplicity $m=1$. Here we drop the index $\mu$, denoting atomic positions as ${\bm{\mathrm{R}}}_{j}^{W}$. This case is equivalent to atomic positions ${\bm{\mathrm{R}}}_{j}^{W}$ forming a Bravais lattice. Note also that multiplicities $m=1$ occur only for symmorphic space groups [43]. Here, plane waves ${Q}_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W})=\mathcal{q}_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W})$ transform according to the one-dimensional IR

[TABLE]

We saw in Eq. (41) that for wave vectors ${\bm{\mathrm{k}}}$ inside the BZ this becomes the identity representation $\Gamma_{1}$. It is illuminating to rederive this result by writing Eq. (42) as

[TABLE]

It follows from Eq. (3) that $\bm{\mathrm{k}}^{\prime}={\bm{\mathrm{g}}}^{-1}\cdot{\bm{\mathrm{k}}}={\bm{\mathrm{k}}}+{\bm{\mathrm{b}}}_{g}$ with a reciprocal lattice vector ${\bm{\mathrm{b}}}_{g}$. Thus Eq. (42) describing the effect of $g$ in real space is equivalent to

[TABLE]

describing the effect of $g$ in reciprocal space. Hence, for $m=1$ plane waves ${Q}_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W})$ transform under $\mathcal{G}_{{\bm{\mathrm{k}}}}$ in a nontrivial way only if the vector ${\bm{\mathrm{k}}}$ is from the border of the BZ when ${\bm{\mathrm{k}}}^{\prime}={\bm{\mathrm{g}}}^{-1}\cdot{\bm{\mathrm{k}}}$ and ${\bm{\mathrm{k}}}$ can differ by a reciprocal lattice vector ${\bm{\mathrm{b}}}_{g}\neq 0$. Otherwise, ${\bm{\mathrm{k}}}^{\prime}={\bm{\mathrm{k}}}$ implies that ${Q}_{\bm{\mathrm{k}}}^{W}({\bm{\mathrm{R}}}_{j}^{W})$ transforms according to $\Gamma_{1}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$.

If a space group has Wyckoff positions $W\neq W^{\prime}$ each with multiplicity $m=1$, we can compare the IRs of the plane waves at the positions ${\bm{\mathrm{R}}}_{j}^{W}$ and ${\bm{\mathrm{R}}}_{j}^{W^{\prime}}$. We have

[TABLE]

Hence

[TABLE]

For $W\neq W^{\prime}$, the vector ${\bm{\mathrm{R}}}_{j}^{W}-{\bm{\mathrm{R}}}_{j}^{W^{\prime}}$ is not equal to a lattice vector, so that for nonzero ${\bm{\mathrm{b}}}_{g}$ (i.e., for ${\bm{\mathrm{k}}}$ on the boundary of the Brillouin zone) we generally have

[TABLE]

This implies that a nontrivial change of the coordinate system which requires a relabeling of the Wyckoff letters associated with atomic positions changes the IRs of the plane waves at these positions. This relabeling of IR assignments will be discussed in Sec. II.6.2.

As a simple example for Eq. (47), consider the case where the positions ${\bm{\mathrm{R}}}_{j}^{W}$ are equal to lattice vectors, i.e., one of the positions ${\bm{\mathrm{R}}}_{j}^{W}$ is at the origin of the coordinate system. Hence, Eq. (44) gives $\Gamma^{q}_{{\bm{\mathrm{k}}}W}=\{\mathcal{D}^{q}_{{\bm{\mathrm{k}}}W}(g)=1:g\in\mathcal{G}_{{\bm{\mathrm{k}}}}\}$, i.e., the plane wave ${Q}_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W})$ transforms according to the identity representation. On the other hand, the IR $\Gamma^{q}_{{\bm{\mathrm{k}}}W^{\prime}}$ at a different Wyckoff position $W^{\prime}$ can never be the identity representation.

Examples of Wyckoff positions with multiplicity $m=1$ are the positions occupied by the Mo atoms and the center of the hexagon in monolayer MoS2. The plane wave at these positions transforms according to different IRs. In a certain coordinate system where one of these two inequivalent positions are located at the origin, the plane waves at that Wyckoff position transform as the identity representation, while the plane waves at the other Wyckoff position transform according to a different IR. Another example is given by the inequivalent IRs of the plane waves at the positions of C atoms in the central layer of graphene with odd number of layers such as the trilayer graphene discussed in Sec. IV.1.3.

II.6.2 Wyckoff positions with multiplicity $m>1$

For Wyckoff positions with multiplicity $m>1$ the simple analysis based on Eq. (3) is not valid as it does not keep track of how positions ${\bm{\mathrm{R}}}_{j}^{W\mu}$ are mapped onto each other by symmetry operations $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$. Instead we need to use the plane-wave spinors (31). For multiplicity $m>1$, the representation $\Gamma^{q}_{{\bm{\mathrm{k}}}W}$ characterizing the plane-wave spinors ${Q}_{\bm{\mathrm{k}}}({\bm{\mathcal{R}}}_{j}^{W})$ is generally reducible. Using the projection operator (122), we can construct linear combinations $\mathcal{Q}_{\bm{\mathrm{k}}}^{I\beta}({\bm{\mathcal{R}}}_{j}^{W})$ of the base spinors $\mathcal{q}_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ transforming like component $\beta$ of the IR $\Gamma_{I}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ contained in $\Gamma^{q}_{{\bm{\mathrm{k}}}W}$

[TABLE]

where we used Eq. (32). This yields

[TABLE]

with expansion coefficients

[TABLE]

that completely characterize each symmetrized spinor $\mathcal{Q}_{\bm{\mathrm{k}}}^{I\beta}({\bm{\mathcal{R}}}_{j}^{W})$.

Upon translation by a lattice vector ${\bm{\mathrm{a}}}$, the plane waves $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ acquire a phase $\exp(i{\bm{\mathrm{k}}}\cdot{\bm{\mathrm{a}}})$

[TABLE]

This implies that $\mathcal{q}_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ represents (for each $\mu$) a discrete Bloch function for wave vector ${\bm{\mathrm{k}}}$. Similarly, linear combinations of these base spinors including the spinors ${Q}_{\bm{\mathrm{k}}}({\bm{\mathcal{R}}}_{j}^{W})$ and $\mathcal{Q}_{\bm{\mathrm{k}}}^{I\beta}({\bm{\mathcal{R}}}_{j}^{W})$ are thus discrete Bloch functions for wave vector ${\bm{\mathrm{k}}}$. The expansion coefficients $u_{{\bm{\mathrm{k}}}W\mu}^{I\beta}$ take the role of lattice-periodic functions for these discrete Bloch functions.

The projection (48) is valid for all wave vectors ${\bm{\mathrm{k}}}$ in the Brillouin zone (though trivial for $m=1$ when ${\bm{\mathrm{k}}}$ is inside the Brillouin zone, as noted above). In general, the projectors $\Pi_{\bm{\mathrm{k}}}^{I\beta}$ decompose a plane-wave spinor ${Q}_{\bm{\mathrm{k}}}({\bm{\mathcal{R}}}_{j}^{W})$ into multiple Bloch functions $\mathcal{Q}_{\bm{\mathrm{k}}}^{I\beta}({\bm{\mathcal{R}}}_{j}^{W})$ corresponding to different IRs $\Gamma_{I}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$. Yet we often have sets of special positions ${\bm{\mathcal{R}}}_{Wj}^{I{\bm{\mathrm{k}}}}=\{{\bm{\mathrm{R}}}_{W\mu j}^{I{\bm{\mathrm{k}}}}:\mu=1,\ldots,m\}$ within the unit cell where only the Bloch functions $\mathcal{Q}_{\bm{\mathrm{k}}}^{I\beta}({\bm{\mathcal{R}}}_{j}^{W})$ for one IR $\Gamma_{I}$ are nonzero, but all other projections vanish. This greatly simplifies further discussion of TB Bloch functions at positions ${\bm{\mathcal{R}}}_{Wj}^{I{\bm{\mathrm{k}}}}$. The positions ${\bm{\mathcal{R}}}_{Wj}^{I{\bm{\mathrm{k}}}}$ are characterized by two different groups, the group describing the site symmetry [43] denoted as $\mathcal{G}_{W}$ and the group of the wave vector $\mathcal{G}_{{\bm{\mathrm{k}}}}$. Often the positions ${\bm{\mathcal{R}}}_{Wj}^{I{\bm{\mathrm{k}}}}$ with nontrivial $\mathcal{G}_{{\bm{\mathrm{k}}}}$ coincide with Wyckoff positions with a nontrivial $\mathcal{G}_{W}$. For Wyckoff positions with multiplicity $m=1$, we always have $\mathcal{G}_{W}=\mathcal{G}_{0}$, so that $\mathcal{G}_{{\bm{\mathrm{k}}}}\subseteq\mathcal{G}_{W}$. However, for $m>1$ we will find below that in general there is no simple relation between the group $\mathcal{G}_{W}$ characterizing such special positions ${\bm{\mathcal{R}}}_{Wj}^{I{\bm{\mathrm{k}}}}$ and the group of the wave vector $\mathcal{G}_{{\bm{\mathrm{k}}}}$ for which this happens 333In SLG, the carbon atoms with multiplicity $m=2$ are characterized by the site symmetry group $\mathcal{G}_{W}=D_{3h}$, see Table 1. Here, likewise, the group of the wave vector at the ${\bm{\mathrm{K}}}$ point is $\mathcal{G}_{{\bm{\mathrm{K}}}}=D_{3h}$ and symmetrized plane waves at ${\bm{\mathrm{K}}}$ transform according to the two-dimensional IR $\Gamma_{6}$ of $D_{3h}$, see Table 6. On the other hand, in BLG, the carbon atoms at positions ${\bm{\mathrm{R}}}_{j}^{d1},{\bm{\mathrm{R}}}_{j}^{d2}$ have site symmetry $\mathcal{G}_{W}=C_{3v}$, whereas the symmetrized plane waves at ${\bm{\mathrm{K}}}$ for these positions transform according to the two-dimensional IR $\Gamma_{3}$ of $D_{3}$. Neither of the groups $C_{3v}$ and $D_{3}$ can be viewed as a subgroup of the other one..

The symmetrized plane waves $\mathcal{Q}_{\bm{\mathrm{k}}}^{I\beta}({\bm{\mathcal{R}}}_{j}^{W})$ including the positions ${\bm{\mathcal{R}}}_{Wj}^{I{\bm{\mathrm{k}}}}$ are universal features of each space group, independent of the “atomistic realization” of a space group in different crystal structures (e.g., the number and positions of atoms in a unit cell). They apply both to spinless models and models that include the spin degree of freedom. We note that the symmetrized plane waves $\mathcal{Q}_{\bm{\mathrm{k}}}^{I\beta}({\bm{\mathcal{R}}}_{j}^{W})$ introduced here in the context of the TB approximation for Wyckoff positions ${\bm{\mathcal{R}}}_{j}^{W}$ are conceptually different from the symmetrized plane waves discussed previously in the context of the nearly-free electron approximation, see, e.g., Refs. [38, 50, 39].

II.7 Rearrangement of IRs of Bloch states under a change of the coordinate system

The Bloch states in certain crystal structures are characterized by IRs of the group $\mathcal{G}_{{\bm{\mathrm{k}}}}$ that depend in a nontrivial way on the location of the origin or the orientation of the coordinate system relative to the position of the atoms [20, 21, 22, 23]. Cornwell [22] has given a general discussion of the origin dependence of the symmetry labeling of electron states in such systems. Here we review and extend these findings, focusing on symmorphic space groups and adopting a notation matching other parts of this study. We show that for different choices of the origin we get a rearrangement of the IRs of Bloch states. We exploit this rearrangement lemma when discussing band symmetries for specific materials further below.

We consider a crystal structure with space group $\mathcal{G}$. For the coordinate system ${\bm{\mathrm{r}}}=(x,y,z)$, the lattice-periodic (single-electron) Hamiltonian is $H({\bm{\mathrm{r}}})$. The eigenfunctions of $H({\bm{\mathrm{r}}})$ are Bloch functions $\Psi_{n{\bm{\mathrm{k}}}}^{I\beta}({\bm{\mathrm{r}}})$, obeying the eigenvalue equation

[TABLE]

with energy $E_{n}({\bm{\mathrm{k}}})$. For a given wave vector ${\bm{\mathrm{k}}}$, the index $I$ denotes the IR $\Gamma_{I}$ of the point group $\mathcal{G}_{{\bm{\mathrm{k}}}}$ of the wave vector, to be discussed in more detail below. The eigenfunction $\Psi_{n{\bm{\mathrm{k}}}}^{I\beta}({\bm{\mathrm{r}}})$ transforms according to the $\beta$th component of the IR $\Gamma_{I}$. For brevity, we drop in this section the band index $n$.

We denote coordinate transformations using the Seitz notation as $\{g\mskip 1.19998mu|\mskip 1.19998mu{\bm{\mathrm{\tau}}}\}$, where $g$ is a (proper or improper) rotation that is followed by a translation ${\bm{\mathrm{\tau}}}$. We seek to identify a pure translation $T\equiv\{\openone\mskip 1.19998mu|\mskip 1.19998mu{\bm{\mathrm{\tau}}}\}$ of the coordinate system ${\bm{\mathrm{r}}}$, where ${\bm{\mathrm{\tau}}}$ equals a fraction of a lattice vector such that for the shifted, primed coordinated system ${\bm{\mathrm{r}}}^{\prime}=(x^{\prime},y^{\prime},z^{\prime})$ the crystal structure has the same space group symmetry $\mathcal{G}$ as for the unprimed coordinate system ${\bm{\mathrm{r}}}=(x,y,z)$. The translation $T$ transforms the Bloch function $\Psi_{{\bm{\mathrm{k}}}}^{I\beta}({\bm{\mathrm{r}}})$ into

[TABLE]

The index $I^{\prime}\neq I$ will be justified below. As $T=\{\openone\mskip 1.19998mu|\mskip 1.19998mu{\bm{\mathrm{\tau}}}\}$ commutes with primitive translations $\{\openone\mskip 1.19998mu|\mskip 1.19998mu{\bm{\mathrm{a}}}\}$ we get

[TABLE]

so that the transformed Bloch function $\Psi_{{\bm{\mathrm{k}}}}^{I^{\prime}\beta}({\bm{\mathrm{r}}}^{\prime})$ has, indeed, the same wave vector ${\bm{\mathrm{k}}}$ as $\Psi_{{\bm{\mathrm{k}}}}^{I\beta}({\bm{\mathrm{r}}})$. Furthermore, the Hamiltonian in the primed coordinate system becomes

[TABLE]

Thus

[TABLE]

so that the transformed function $\Psi_{{\bm{\mathrm{k}}}}^{I^{\prime}\beta}({\bm{\mathrm{r}}}^{\prime})$ is an eigenfunction of $H({\bm{\mathrm{r}}}^{\prime})$ with the same eigenvalue $E({\bm{\mathrm{k}}})$ as $\Psi_{{\bm{\mathrm{k}}}}^{I\beta}({\bm{\mathrm{r}}})$.

Given the space group $\mathcal{G}$ of the crystal, the invariance of $H({\bm{\mathrm{r}}})$ under the symmetry operations $\mathcal{T}\equiv\{g\mskip 1.19998mu|\mskip 1.19998mu{\bm{\mathrm{a}}}\}\in\mathcal{G}$ reads

[TABLE]

For the Hamiltonian $H({\bm{\mathrm{r}}}^{\prime})$ in the primed coordinate system we get

[TABLE]

Hence the primed Hamiltonian $H({\bm{\mathrm{r}}}^{\prime})$ obeys an invariance condition analogous to Eq. (57) if $(T^{-1}\,\mathcal{T}\,T)$ commutes with $H({\bm{\mathrm{r}}})$

[TABLE]

This requires in turn, given Eq. (57), that for all $\mathcal{T}\in\mathcal{G}$

[TABLE]

is an element of the space group $\mathcal{G}$, so that $g\,{\bm{\mathrm{\tau}}}-{\bm{\mathrm{\tau}}}$ must be equal to a lattice vector ${\bm{\mathrm{a}}}^{\prime}$ of the crystal

[TABLE]

where $\mathcal{G}_{0}$ is the point group corresponding to the space group $\mathcal{G}$. A nontrivial solution to this problem is a vector ${\bm{\mathrm{\tau}}}$ that is not equal to a lattice vector ${\bm{\mathrm{a}}}$. Equation (61) defines the allowed shifts ${\bm{\mathrm{\tau}}}$ (up to a lattice vector) that provide alternative descriptions of a crystal structure with space group $\mathcal{G}$.

We obtain nontrivial solutions to Eq. (61), for example, if a crystal consists of atoms at Wyckoff positions $W=A,B,\ldots$ each with multiplicity $m=1$ 444In Cornwell’s notation [22], different Wyckoff positions each with multiplicity $m=1$ are positions each forming a Bravais lattice.. We denote these positions in the unit cell as ${\bm{\mathrm{t}}}_{\mathrm{A}}$, ${\bm{\mathrm{t}}}_{\mathrm{B}}$, ${\bm{\mathrm{t}}}_{\mathrm{C}}$, $\dots$, respectively. By definition of the space group $\mathcal{G}$, these positions ${\bm{\mathrm{t}}}_{W}$ obey the condition

[TABLE]

with lattice vectors ${\bm{\mathrm{a}}}^{\prime}$. Hence, any linear combination of these position vectors $\{{\bm{\mathrm{t}}}_{W}\}$ with integer prefactors (e.g., ${\bm{\mathrm{\tau}}}={\bm{\mathrm{t}}}_{W}-{\bm{\mathrm{t}}}_{W^{\prime}}$ with $W\neq W^{\prime}$) yields a translation ${\bm{\mathrm{\tau}}}$ consistent with Eq. (61).

In the unprimed coordinate system the eigenfunctions $\Psi_{{\bm{\mathrm{k}}}}^{I\beta}({\bm{\mathrm{r}}})$ transform according to the $\beta$th component of an IR $\Gamma_{I}$ of the point group $\mathcal{G}_{{\bm{\mathrm{k}}}}$ of the wave vector ${\bm{\mathrm{k}}}$

[TABLE]

with representation matrices $\mathcal{D}_{I}(g)_{\beta^{\prime}\beta}$. We can evaluate the action of a symmetry operation $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$ on primed Bloch functions $\Psi_{{\bm{\mathrm{k}}}}^{I^{\prime}\beta}({\bm{\mathrm{r}}}^{\prime})$ as follows

[TABLE]

where the primed representation matrices $\mathcal{D}_{I^{\prime}}(g)_{\beta^{\prime}\beta}$ become

[TABLE]

We denote the phase factors in Eq. (65) by

[TABLE]

Using Eq. (3), this becomes

[TABLE]

The phase $\mathcal{D}_{\bm{\mathrm{\tau}}}^{\bm{\mathrm{k}}}(g)$ is therefore nontrivial when ${\bm{\mathrm{b}}}_{g}\neq 0$, which can only happen at the border of the Brillouin zone.

We show in the next paragraph that $\Gamma_{\bm{\mathrm{\tau}}}\equiv\{\mathcal{D}_{\bm{\mathrm{\tau}}}^{\bm{\mathrm{k}}}(g):g\in\mathcal{G}_{{\bm{\mathrm{k}}}}\}$ defines a one-dimensional IR of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ (for every wave vector ${\bm{\mathrm{k}}}$ in the Brillouin zone). Therefore, the IR $\Gamma_{I^{\prime}}$ of a Bloch function in the primed coordinate system is given by $\Gamma_{\bm{\mathrm{\tau}}}$ times the IR $\Gamma_{I}$ of the Bloch function in the unprimed coordinate system, so that Eq. (65) becomes

[TABLE]

The rearrangement lemma for IRs discussed in Appendix B applied to Eq. (68) shows that, unless we have the trivial case that $\Gamma_{\bm{\mathrm{\tau}}}$ is the identity representation, each IR $\Gamma_{I}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ in the unprimed coordinate system is mapped on an IR $\Gamma_{I^{\prime}}\neq\Gamma_{I}$ in the primed coordinate system. Hence we call Eq. (68) the rearrangement lemma for the IRs of Bloch states and $\Gamma_{\bm{\mathrm{\tau}}}$ the rearrangement representation (RAR) for the coordinate shift ${\bm{\mathrm{\tau}}}$. Examples for this rearrangement of IRs of Bloch states will be given below when we study the symmetries of the Bloch functions in monolayer MoS2 (Sec. III.2) and trilayer graphene (Sec. IV.1.3). It follows from Eq. (67) that only at the border of the Brillouin zone the IR labeling of Bloch states can depend on the origin of the coordinate system [22]. Also, Eq. (67) implies $\Gamma_{-{\bm{\mathrm{\tau}}}}=\Gamma_{\bm{\mathrm{\tau}}}^{\ast}$. Generally, the shift ${\bm{\mathrm{\tau}}}$ is defined up to a lattice vector ${\bm{\mathrm{a}}}$. It follows immediately from Eq. (67) that ${\bm{\mathrm{\tau}}}$ and $\tilde{{\bm{\mathrm{\tau}}}}\equiv{\bm{\mathrm{\tau}}}+{\bm{\mathrm{a}}}$ define the same RAR $\Gamma_{\bm{\mathrm{\tau}}}$.

To show that $\Gamma_{\bm{\mathrm{\tau}}}=\{\mathcal{D}_{\bm{\mathrm{\tau}}}^{\bm{\mathrm{k}}}(g):g\in\mathcal{G}_{{\bm{\mathrm{k}}}}\}$ defines a one-dimensional IR of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ we consider two group elements $g_{i}\in\mathcal{G}_{{\bm{\mathrm{k}}}}$ ($i=1,2$). According to Eq. (61), the transformations $g_{i}\,{\bm{\mathrm{\tau}}}$ differ from ${\bm{\mathrm{\tau}}}$ by lattice vectors ${\bm{\mathrm{a}}}_{i}$

[TABLE]

so that

[TABLE]

and

[TABLE]

We can also write $g_{1}g_{2}\,{\bm{\mathrm{\tau}}}$ as

[TABLE]

Hence

[TABLE]

We get similar to Eq. (2) and using Eq. (3)

[TABLE]

This gives us finally

[TABLE]

so that, indeed, $\Gamma_{\bm{\mathrm{\tau}}}=\{\mathcal{D}_{\bm{\mathrm{\tau}}}^{\bm{\mathrm{k}}}(g):g\in\mathcal{G}_{{\bm{\mathrm{k}}}}\}$ is a one-dimensional IR of the group $\mathcal{G}_{{\bm{\mathrm{k}}}}$.

The above argument [22] is based on the full Bloch functions $\Psi_{\bm{\mathrm{k}}}({\bm{\mathrm{r}}})$. We showed in Eq. (16) that in a TB description, these Bloch functions can be factorized as $\Phi_{\nu{\bm{\mathrm{k}}}}^{W\mu}({\bm{\mathrm{r}}})=q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})\,\phi_{\nu}^{W}(\delta{\bm{\mathrm{r}}})$, so that the symmetry of the plane waves $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ can be discussed separately from the symmetry of the atomic orbitals $\phi_{\nu}^{W}(\delta{\bm{\mathrm{r}}})$. The orbitals $\phi_{\nu}^{W}(\delta{\bm{\mathrm{r}}})$ only depend on $\delta{\bm{\mathrm{r}}}$ but not on the actual positions ${\bm{\mathrm{R}}}_{j}^{W\mu}$. Hence it follows immediately that the symmetry of the orbitals $\phi_{\nu}^{W}(\delta{\bm{\mathrm{r}}})$ is independent of the coordinate system used, see also Eq. (26). Only the representation $\Gamma^{q}_{{\bm{\mathrm{k}}}W}$ of the plane waves $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ depends, in general, on the coordinate system. A translation of the coordinate system by ${\bm{\mathrm{\tau}}}$ maps the positions ${\bm{\mathrm{R}}}_{j}^{W\mu}$ onto the atomic position ${\bm{\mathrm{R}}}_{j}^{W^{\prime}\mu}={\bm{\mathrm{R}}}_{j}^{W\mu}-{\bm{\mathrm{\tau}}}$. The new coordinate system is valid if and only if $g\,{\bm{\mathrm{\tau}}}-{\bm{\mathrm{\tau}}}$ are lattice vectors for all $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$ [see Eq. (61)], so that the plane waves transform as $\Gamma^{q}_{{\bm{\mathrm{k}}}W^{\prime}}$. We have

[TABLE]

so that for the transformed Wyckoff letter $W^{\prime}$, Eq. (33) becomes

[TABLE]

with $\mathcal{D}_{\bm{\mathrm{\tau}}}^{\bm{\mathrm{k}}}(g)$ given by Eq. (66). This gives us the rearrangement lemma for the representations of plane waves

[TABLE]

Hence we confirm that the nontrivial case $\Gamma^{q}_{{\bm{\mathrm{k}}}W^{\prime}}\neq\Gamma^{q}_{{\bm{\mathrm{k}}}W}$ requires that $\Gamma_{\bm{\mathrm{\tau}}}=\{\mathcal{D}_{\bm{\mathrm{\tau}}}^{\bm{\mathrm{k}}}(g):g\in\mathcal{G}_{{\bm{\mathrm{k}}}}\}$ is not the identity representation. As mentioned above, the symmetry of plane waves $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{W\mu})$ is a universal problem for each space group $\mathcal{G}$, independent of the detailed realization of a crystal structure. This holds, in particular, if the atoms are located at positions ${\bm{\mathrm{R}}}_{W\mu j}^{I{\bm{\mathrm{k}}}}$ where the $\mathcal{q}_{\bm{\mathrm{k}}}^{I\beta}({\bm{\mathrm{R}}}_{W\mu j}^{I{\bm{\mathrm{k}}}})$ transforms according to only one IR $\Gamma_{I}$. Hence it is possible to discuss the rearrangement lemma for the IRs of Bloch states independent of a particular crystal structure, but it depends only on the space group $\mathcal{G}$. Among all 230 space groups, 159 contain Wyckoff sites with origin-dependent site symmetries [52]. Though a necessary criterion, it is however not a sufficient criterion for a rearrangement of the IRs of Bloch states under a change of the coordinate system [22].

II.8 Effect of time reversal

In the absence of an external magnetic field, the eigenfunctions of the TB Hamiltonian obey time-reversal symmetry $\Theta$. This implies that if an eigenfunction $\Psi_{{\bm{\mathrm{k}}}}^{I\beta}({\bm{\mathrm{r}}})$ with energy $E({\bm{\mathrm{k}}})$ transforms according to the $\beta$th component of the IR $\Gamma_{I}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$, the time-reversed wave function $\Theta\,\Psi_{{\bm{\mathrm{k}}}}^{I\beta}({\bm{\mathrm{r}}})=\Psi_{{\bm{\mathrm{k}}}}^{I\beta}({\bm{\mathrm{r}}})^{\ast}=\Psi_{-{\bm{\mathrm{k}}}}^{\tilde{I},\beta^{\prime}}({\bm{\mathrm{r}}})$ [which is likewise an eigenfunction of the Hamiltonian with energy $E(-{\bm{\mathrm{k}}})=E({\bm{\mathrm{k}}})$] transforms according to the $\beta^{\prime}$th component of the complex conjugate IR $\Gamma_{I}^{\ast}$ of $\mathcal{G}_{-{\bm{\mathrm{k}}}}=\mathcal{G}_{{\bm{\mathrm{k}}}}$. Therefore, if the eigenfunctions of the TB Hamiltonian for some energy $E({\bm{\mathrm{k}}})$ contain atomic orbitals transforming according to an IR $\Gamma_{J}^{\phi}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ and symmetry-adapted plane waves transforming according to $\Gamma_{J^{\prime}}^{q}$, the eigenfunctions for wave vector $-{\bm{\mathrm{k}}}$ with energy $E(-{\bm{\mathrm{k}}})=E({\bm{\mathrm{k}}})$ contain atomic orbitals transforming according to the complex conjugate IR $\Gamma_{J}^{\phi\ast}$ and plane waves transforming according to $\Gamma_{J^{\prime}}^{q\ast}$. Here, the symmetry-adapted atomic orbitals at $-{\bm{\mathrm{k}}}$ are the complex conjugates of the corresponding atomic orbitals at ${\bm{\mathrm{k}}}$. We obtain the symmetry-adapted plane waves at $-{\bm{\mathrm{k}}}$ from the corresponding plane waves at ${\bm{\mathrm{k}}}$ by replacing ${\bm{\mathrm{k}}}\rightarrow-{\bm{\mathrm{k}}}$.

Degeneracies of Bloch states due to time-reversal symmetry are discussed in Refs. [53, 38, 3]. In general, three cases must be distinguished 555The classification of representations under time reversal adopted here from Ref. [3] is the most convenient for physical applications, but it differs from that customary in courses on group theory, where real representations are assigned to case (a) and complex inequivalent and equivalent representations to cases (b) and (c). The two classifications agree for single-group representations, but cases (a) and (c) must be interchanged for double-group representations.. In case (a), eigenfunctions $\Psi$ and $\Theta\,\Psi$ of the crystal Hamiltonian $H({\bm{\mathrm{r}}})$ are linearly dependent. In case (b), eigenfunctions $\Psi$ and $\Theta\,\Psi$ of $H({\bm{\mathrm{r}}})$ are linearly independent and transform according to inequivalent representations $\Gamma_{I}$ and $\Gamma_{I}^{\ast}$, i.e., $\chi_{I}(g)\neq\chi_{I}^{\ast}(g)$ for some $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$. Finally, in case (c), eigenfunctions $\Psi$ and $\Theta\,\Psi$ of $H({\bm{\mathrm{r}}})$ are linearly independent and transform according to equivalent representations $\Gamma_{I}$ and $\Gamma_{I}^{\ast}$, i.e., $\chi_{I}(g)=\chi_{I}^{\ast}(g)$ for all $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$. In cases (b) and (c) invariance under time reversal causes additional degeneracy. We have for symmorphic space groups [53, 38, 3]

[TABLE]

where $f$ is the number of points of the star of ${\bm{\mathrm{k}}}$, $h$ is the order of the crystallographic point group $\mathcal{G}_{0}$ of the crystal, $\chi_{I}(g)$ are the characters of the IR $\Gamma_{I}$ of the point group $\mathcal{G}_{{\bm{\mathrm{k}}}}$ of the wave vector ${\bm{\mathrm{k}}}$, and $g{\bm{\mathrm{k}}}+{\bm{\mathrm{k}}}$ may be zero or a reciprocal lattice vector ${\bm{\mathrm{b}}}$. We have $\Theta^{2}=+1$ for single-group representations and $\Theta^{2}=-1$ for double-group representations. The criterion (79) applies to the IRs $\Gamma_{I}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ independent of the origin of the coordinate system. If a crystal structure permits a change of the coordinate system characterized by a vector ${\bm{\mathrm{\tau}}}$ with RAR $\Gamma_{\bm{\mathrm{\tau}}}$, a Bloch function transforming according to the IR $\Gamma_{I}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ in the old coordinate system transforms according to $\Gamma_{I^{\prime}}=\Gamma_{\bm{\mathrm{\tau}}}\times\Gamma_{I}$ in the new coordinate, see Eq. (68). Therefore, the IRs $\Gamma_{I}$ and $\Gamma_{I^{\prime}}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ must fall into the same category according to Eq. (79).

In a more detailed analysis [55, 38, 3], for each of the cases (a), (b), and (c) three possibilities must be distinguished: (1) the points ${\bm{\mathrm{k}}}$ and $-{\bm{\mathrm{k}}}$ are equivalent, i.e., ${\bm{\mathrm{k}}}=-{\bm{\mathrm{k}}}+{\bm{\mathrm{b}}}$; (2) ${\bm{\mathrm{k}}}$ is not equivalent to $-{\bm{\mathrm{k}}}$, but the space group contains an element $R$ which maps ${\bm{\mathrm{k}}}$ onto $-{\bm{\mathrm{k}}}$

[TABLE]

(3) the points ${\bm{\mathrm{k}}}$ and $-{\bm{\mathrm{k}}}$ are in different stars. For the systems discussed below, case (1) applies to the ${\bm{\mathrm{\Gamma}}}$ and ${\bm{\mathrm{M}}}$ points of the BZ, whereas case (2) applies to the ${\bm{\mathrm{K}}}$ points.

III Band Symmetries in MoS2

Having derived a systematic theory for the symmetry of TB Bloch functions, we now apply this theory to several quasi-2D materials. Our main focus is on monolayer MoS2. For comparison, we also discuss single-layer (SLG), bilayer (BLG), and trilayer (TLG) graphene in the next section.

III.1 Crystal structure of MoS2

The crystal structure of single-layer MoS2 is shown in Fig. 1. It is characterized by the point group $D_{3h}$ and space group $P\bar{6}m2$ (# 187). Three Wyckoff positions have multiplicity $m=1$, the positions of the Mo atom, the midpoint between a pair of top and bottom S atoms, and the center of the hexagon. Hence, as discussed in Sec. II.7, three choices $\alpha=a,b,c$ emerge for the origin of the coordinate system: ($a$) origin at the center of the hexagon [Fig. 1(a)], ($b$) origin at a Mo atom [Fig. 1(b)], and ($c$) origin at the midpoint between a top and bottom S atom [Fig. 1(c)]. In either case, the Mo atoms have site symmetry $D_{3h}$ and Wyckoff multiplicity $m=1$. Yet the Wyckoff letters for these positions listed in Table 1 depend on the coordinate system [52]. For coordinate system ($a$), the Mo and S atoms have Wyckoff letters $e$ and $h$, respectively, whereas these letters become $a$ and $i$ in coordinate system ($b$), and $c$ and $g$ in coordinate system ($c$). The S atoms have site symmetry $C_{3v}$ and multiplicity $m=2$. The positions of Mo and S atoms in unit cell $j$ are denoted by ${\bm{\mathrm{R}}}_{j}^{\textrm{Mo}(\alpha)}$ and ${\bm{\mathrm{R}}}_{j}^{\textrm{S}\mu(\alpha)}$, respectively. For the S atom, the top (bottom) atoms are labeled $\mu=1$ ($\mu=2$). There are yet other coordinate systems that can be used for MoS2. For example, Ref. [31] used a coordinate system that differs from coordinate system ($a$) by a reflection about the $xz$ plane. Here, we do not consider these additional coordinate systems [23].

The Brillouin zone for single-layer MoS2 is shown in Fig. 1(e). In the following, we will focus on the ${\bm{\mathrm{\Gamma}}}$ point ${\bm{\mathrm{k}}}=0$, the ${\bm{\mathrm{K}}}$, and the ${\bm{\mathrm{M}}}$ points. The star of the ${\bm{\mathrm{K}}}$ point includes two inequivalent wave vectors denoted ${\bm{\mathrm{K}}}$ and ${\bm{\mathrm{K}}}^{\prime}$. The star of the ${\bm{\mathrm{M}}}$ point includes three inequivalent wave vectors denoted ${\bm{\mathrm{M}}}_{1}$, ${\bm{\mathrm{M}}}_{2}$, and ${\bm{\mathrm{M}}}_{3}$.

The point group of single-layer MoS2 is $D_{3h}$. This is also the point group of the wave vector at the ${\bm{\mathrm{\Gamma}}}$ point. It contains a $120^{\circ}$ counterclockwise rotation $C_{3}$ about the $z$ axis. The reflection plane of $\sigma_{h}$ is the $xy$ plane and $S_{3}=\sigma_{h}C_{3}$. The rotation axes of the three twofold rotations $C_{2}^{\prime(i)}$ are the axes $i=1,2,3$ shown in Fig. 1. These axes are also shown as dashed lines in Fig. 2(a). The reflection plane of $\sigma_{v}^{(i)}$ is the plane passing through the axis $i$ and the $z$ axis. The characters of $D_{3h}$ are listed in Table 25 (Appendix C). We label the IRs of the crystallographic point groups following Koster et al. [56].

At the ${\bm{\mathrm{K}}}$ point [Fig. 2(b)], the point group of the wave vector becomes $\mathcal{G}_{{\bm{\mathrm{K}}}}=C_{3h}$ whose characters are listed in Table 27 (Appendix C). Finally, at the inequivalent points ${\bm{\mathrm{M}}}_{i}$ ($i=1,2,3$), the group of the wave vector is $\mathcal{G}_{{\bm{\mathrm{M}}}}=C_{2v}$, the character table of which is reproduced in Table 28 (Appendix C). This group contains the twofold rotation $C_{2}$, the axis of which is indicated as dashed line in Figs. 2(c)–2(e), the reflection $\sigma_{v}$ about the $xy$ plane, and the reflection $\sigma_{v}^{\prime}$ for which the reflection plane includes the dashed line and the $z$ axis.

The primitive lattice vectors ${\bm{\mathrm{a}}}_{1}$ and ${\bm{\mathrm{a}}}_{2}$ are

[TABLE]

where $a$ is the lattice constant. Ignoring for brevity the $z$ component, the positions ${\bm{\mathrm{t}}}_{\mathrm{Mo}}^{\alpha}$ of Mo and ${\bm{\mathrm{t}}}_{\mathrm{S}}^{\alpha}$ of S in the unit cell are

[TABLE]

where the superscript $\alpha=a,b,c$ denotes the coordinate system. The primitive vectors ${\bm{\mathrm{b}}}_{1}$ and ${\bm{\mathrm{b}}}_{2}$ of the reciprocal lattice are

[TABLE]

The two inequivalent corner points of the Brillouin zone are

[TABLE]

and the ${\bm{\mathrm{M}}}$ points are

[TABLE]

see Fig. 1(e).

III.2 Rearrangement of IRs of Bloch states in MoS2

The crystal structure of monolayer MoS2 can be described by three different coordinate systems $\alpha=a,b,c$ shown in Fig. 1. This provides an example for the rearrangement of the IRs of Bloch states discussed in general terms in Sec. II.7. The coordinate systems $\alpha$ and $\beta$ are related via a translation ${\bm{\mathrm{\tau}}}_{\alpha\beta}$. For the shifts $a\rightarrow b$, $b\rightarrow c$, and $c\rightarrow a$, the translation vectors (apart from a lattice vector) are given by

[TABLE]

and ${\bm{\mathrm{\tau}}}_{\beta\alpha}=-{\bm{\mathrm{\tau}}}_{\alpha\beta}$.

For wave vectors ${\bm{\mathrm{k}}}$ inside the BZ such as the ${\bm{\mathrm{\Gamma}}}$ point as well as for the ${\bm{\mathrm{M}}}$ points, the RARs $\Gamma_{\alpha\beta}^{\bm{\mathrm{\Gamma}}}$ and $\Gamma_{\alpha\beta}^{\bm{\mathrm{M}}}$ are given by Eq. (67) with ${\bm{\mathrm{b}}}_{g}=0$ for all $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$. These RARs are, therefore, given by the identity representation $\Gamma_{1}$, i.e., at both the ${\bm{\mathrm{\Gamma}}}$ and ${\bm{\mathrm{M}}}$ points the labeling of Bloch states is independent of the coordinate system. However, a shift of the coordinate system rearranges the IRs at the ${\bm{\mathrm{K}}}$ point. Using Eq. (67b) at the ${\bm{\mathrm{K}}}$ point, where the group of the wave vector is $C_{3h}$, we get

[TABLE]

For the shifts $a\rightarrow b$, $b\rightarrow c$, and $c\rightarrow a$, we thus have using Eq. (67)

[TABLE]

with $\omega\equiv\exp(i\pi/6)$. This implies $\Gamma_{ab}^{\bm{\mathrm{K}}}=\Gamma_{bc}^{\bm{\mathrm{K}}}=\Gamma_{ca}^{\bm{\mathrm{K}}}=\Gamma_{3}$. Since ${\bm{\mathrm{\tau}}}_{\beta\alpha}=-{\bm{\mathrm{\tau}}}_{\alpha\beta}$, we have $\Gamma_{ba}^{\bm{\mathrm{K}}}=\Gamma_{cb}^{\bm{\mathrm{K}}}=\Gamma_{ac}^{\bm{\mathrm{K}}}=\Gamma_{3}^{\ast}=\Gamma_{2}$. Hence, at the ${\bm{\mathrm{K}}}$ point, for a Bloch state transforming in one coordinate system according to a certain IR, we can multiply this IR with either $\Gamma_{2}$ or $\Gamma_{3}$ to obtain the IR of the same Bloch state in a different coordinate system. The multiplication table for the IRs of $C_{3h}$ is reproduced in Table 30 (Appendix 29). Table 2 shows how the IRs of $C_{3h}$ are rearranged when we go from one of the three coordinate systems $\alpha=a,b,c$ to a different coordinate system $\beta$. At the ${\bm{\mathrm{K}}}^{\prime}$ point, Eq. (66) gives

[TABLE]

III.3 Atomic orbitals $\phi_{\nu}^{W}(\delta\mathbf{r})$ at $\mathbf{k}=\mathbf{\Gamma}$, $\mathbf{K}$, and $\mathbf{M}$

The conduction and valence bands in MoS2 are dominated by Mo $d$ and S $p$ orbitals [24]. At the points ${\bm{\mathrm{\Gamma}}}$, ${\bm{\mathrm{K}}}$, and ${\bm{\mathrm{M}}}$, the groups of the wave vectors $\mathcal{G}_{{\bm{\mathrm{k}}}}$ are $D_{3h}$, $C_{3h}$, and $C_{2v}$, respectively, with character tables reproduced in Tables 25, 27, and 28. We use Eq. (122) to project these functions onto functions transforming according to the IRs of the various groups $\mathcal{G}_{{\bm{\mathrm{k}}}}$. The relevant symmetry operations are defined in Fig. 2 and Table 3 considering both polar vectors ${\bm{\mathrm{P}}}$ such as position ${\bm{\mathrm{r}}}$ and axial vectors ${\bm{\mathrm{A}}}$. The two inequivalent points ${\bm{\mathrm{K}}}$ and ${\bm{\mathrm{K}}}^{\prime}$ are related by a vertical reflection that transforms the component $x$ into $-x$ while keeping the $y$ and $z$ components fixed. Table 4 summarizes the IRs of the $p$ and $d$ orbitals. We note that these results are fully consistent with the full rotation group compatibility tables in Ref. [56].

III.4 Transformation of plane waves $q_{\mathbf{k}}(\mathbf{R}_{j}^{\mathrm{Mo}(\alpha)})$

We now determine the IRs of the plane waves $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{\textrm{Mo}(\alpha)})$ for the coordinate systems $\alpha=a,b,c$ at the ${\bm{\mathrm{\Gamma}}}$, ${\bm{\mathrm{K}}}$, and ${\bm{\mathrm{M}}}$ points of the Brillouin zone. Since the Wyckoff letter corresponding to the positions of the Mo atoms has multiplicity $m=1$, we can use either Eq. (42) or Eq. (44) to determine the phase $\mathcal{D}_{\bm{\mathrm{k}}}^{\mathrm{Mo}(\alpha)}(g)$ acquired by the plane waves $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{\textrm{Mo}(\alpha)})$ under a transformation $g$. We can then derive the IRs of the plane waves using the projection operators (122). The results are summarized in Table 6.

III.5 Transformation of plane waves $q_{\mathbf{k}}(\mathbf{R}_{j}^{\mathrm{S}\mu(\alpha)})$

The S atoms are located at Wyckoff positions with multiplicity $m=2$, so that we represent the plane wave at the positions ${\bm{\mathrm{R}}}_{j}^{\mathrm{S}\mu(\alpha)}$ as a two-component spinor

[TABLE]

We can then use Eq. (33) to determine the phases acquired under symmetry transformations. To obtain the plane wave IRs, it is again advantageous to consider the simplest coordinate system. For the S atoms, this is coordinate system $(c)$ where the origin of the coordinate system is at the midpoint between the S atoms in the top and bottom layer of a unit cell. In this case, the transformation $g$ maps ${\bm{\mathrm{R}}}_{1}^{\mathrm{S}\mu(c)}$ either onto itself or onto ${\bm{\mathrm{R}}}_{1}^{\mathrm{S}\mu^{\prime}(c)}$ with $\mu\neq\mu^{\prime}$, so that Eq. (33) becomes

[TABLE]

for all $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$. We can then determine the IRs of the plane waves for the coordinate systems $(a)$ and $(b)$ by using the respective RAR derived in Sec. III.2. The results are summarized in Table 6.

III.6 IRs of Bloch states in MoS2

The full symmetry-adapted Bloch functions are written as products of symmetrized plane waves and symmetrized atomic orbitals. The five symmetry-adapted $d$ orbitals of the Mo atom times the plane wave $q_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{\mathrm{Mo}(\alpha)})$ and the three symmetry-adapted $p$ orbitals of the S atoms times the plane waves $\mathcal{q}_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{\mathrm{S}1(\alpha)})\pm\mathcal{q}_{\bm{\mathrm{k}}}({\bm{\mathrm{R}}}_{j}^{\mathrm{S}2(\alpha)})$ therefore comprises eleven symmetry-adapted basis functions for MoS2 [30]. The corresponding IRs are listed in Table 7 for the ${\bm{\mathrm{\Gamma}}}$, ${\bm{\mathrm{K}}}$, and ${\bm{\mathrm{M}}}$ points. We list in Tables 8, 9, and 10 the sets of Bloch states transforming according to an IR of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ for the wave vectors ${\bm{\mathrm{k}}}={\bm{\mathrm{\Gamma}}},{\bm{\mathrm{K}}},{\bm{\mathrm{M}}}$. Using these symmetrized Bloch functions, the TB Hamiltonian for a wave vector ${\bm{\mathrm{k}}}$ can be written in a block-diagonal form, where each block refers to the basis functions transforming according to an IR $\Gamma_{I}$ of $\mathcal{G}_{{\bm{\mathrm{k}}}}$ 666The symmetry-adapted bases derived here for MoS2 define unitary transformations for block-diagonalizing a MoS2 TB Hamiltonian that agree with the unitary transformations discussed in Ref. [30].. To classify the additional degeneracy of the Bloch states due to time-reversal symmetry, we evaluate Eq. (79). All IRs of the space group for the stars $\{{\bm{\mathrm{\Gamma}}}\}$, $\{{\bm{\mathrm{K}}}\}$, and $\{{\bm{\mathrm{M}}}\}$ belong to case $(a)$.

IV Band Symmetries in Few-Layer Graphene

We can also apply the general formalism in Sec. II to identify the band symmetries in other quasi-2D materials such as SLG, BLG, and TLG.

IV.1 Crystal structure of few-layer graphene

Like the crystal structure of monolayer MoS2, the crystal structures of SLG, BLG, and TLG belong to the hexagonal crystal system. Therefore, we use the same expressions for the primitive lattice vectors [Eq. (81)] and reciprocal lattice vectors [Eq. (83)]; and we have the same high-symmetry points in the BZ denoted ${\bm{\mathrm{\Gamma}}}$, ${\bm{\mathrm{K}}}$ [Eq. (84)], and ${\bm{\mathrm{M}}}$ [Eq. (92)]. The space groups for few-layer graphene are listed in Table 1. This table also contains the site symmetries of the high-symmetry points for these crystal structures.

IV.1.1 Single-layer graphene

Figure 3 shows the crystal structure of SLG. It is characterized by the point group $D_{6h}$ (space group $P6/mmm$, # 191). The carbon atoms form two distinct Bravais lattices denoted as sublattices $A$ and $B$. The atomic positions denoted by $\{{\bm{\mathrm{R}}}_{j}^{c1},{\bm{\mathrm{R}}}_{j}^{c2}\}$ have site symmetries characterized by the point group $D_{3h}$ and Wyckoff letter $c$. The center of the hexagon, characterized by site symmetry $D_{6h}$ is the only Wyckoff position with multiplicity $m=1$. This point is the origin of the coordinate system for this crystal structure. The positions of the C atoms in the unit cell are given by

[TABLE]

The point group $D_{6h}$ of the crystal, which also characterizes the ${\bm{\mathrm{\Gamma}}}$ point of the BZ, contains twofold, threefold and sixfold rotations $C_{2}$, $C_{3}$ and $C_{6}$ where the $z$ axis is the rotation axis. The rotation axes of the three twofold rotation $C_{2}^{\prime(i)}$ and $C_{2}^{\prime\prime(jk)}$ are the corresponding dashed lines shown in Fig. 4(a). The reflection planes for $\sigma_{v}^{(i)}$ and $\sigma_{d}^{(ij)}$ are perpendicular to the $xy$ plane passing through the corresponding dashed lines. The reflection $\sigma_{h}$ is along the $xy$ plane and $S_{n}=C_{n}\sigma_{h}$. At the ${\bm{\mathrm{K}}}$ point, the group of the wave vector is $\mathcal{G}_{{\bm{\mathrm{K}}}}=D_{3h}$ containing three twofold rotations $C_{2}^{\prime(ij)}$ about the corresponding dashed axes in Fig. 4(b). The reflection plane of $\sigma_{v}^{(ij)}$ is perpendicular to the $xy$ plane along the rotation axis of $C_{2}^{\prime(ij)}$. The threefold rotation axis is the $z$ axis, and $\sigma_{h}$ is a reflection about the $xy$ plane. The group of the wave vector at the ${\bm{\mathrm{M}}}_{i}$ points is $\mathcal{G}_{{\bm{\mathrm{M}}}}=D_{2h}$ containing the symmetry operations $C_{2}$, $C_{2}^{\prime}$, $C_{2}^{\prime\prime}$, $\sigma_{v}$, $\sigma_{v}^{\prime}$, and $\sigma_{v}^{\prime\prime}$ with rotation axes and reflection planes shown in Figs. 4(c)–4(e). The character table for $D_{2h}$ is reproduced in Table 26 (Appendix C).

IV.1.2 Bilayer graphene

The point group $D_{3d}$ (space group $P\bar{3}m1$, # 164) characterizes BLG as shown in Fig. 5. The only Wyckoff position with multiplicity $m=1$ is the midpoint of two C atoms on top of each other. We use this point as the origin of the coordinate system. The atomic positions in BLG are the Wyckoff positions $c$ and $d$, each with multiplicity $m=2$ and site symmetry $C_{3v}$. The two atoms in one Wyckoff letter are labeled $\mu=1$ ($\mu=2$) for the atom in the top (bottom) layer. The layers are arranged in an $AB$ stacking, so that the atomic position ${\bm{\mathrm{R}}}_{j}^{c1}$ is located on top of ${\bm{\mathrm{R}}}_{j}^{c2}$. Ignoring the $z$ component, the positions of the C atoms in the unit cell are

[TABLE]

The threefold proper and sixfold improper rotation axis is the $z$ axis. The axis of the twofold rotation $C_{2}^{\prime(ij)}$ is the corresponding dashed line in Fig. 6(a). The three reflection planes corresponding to $\sigma_{d}^{(i)}$ are perpendicular to the $xy$ plane passing through the corresponding dashed line. The group of the wave vector is $\mathcal{G}_{{\bm{\mathrm{K}}}}=D_{3}$ at the ${\bm{\mathrm{K}}}$ point with threefold rotations about the $z$ axis and three twofold rotations $C_{2}^{\prime(jk)}$ about the corresponding dashed line in Fig. 6(b). At the ${\bm{\mathrm{M}}}_{i}$ points, the group of the wave vector is $\mathcal{G}_{{\bm{\mathrm{M}}}}=C_{2h}$, where $\sigma_{h}$ is perpendicular to the $xy$ plane passing through the corresponding dashed line in Figs. 6(c)–6(e), and the twofold rotation $C_{2}$ is about the corresponding dashed line. The character table for $C_{2h}$ is reproduced in Table 29 (Appendix C).

IV.1.3 Trilayer graphene

Last, Fig. 7 shows the crystal structure of TLG. This system has the same space group $P\bar{6}m2$, # 187 as monolayer MoS2 (point group $D_{3h}$). We designate the Wyckoff positions of the carbon atoms in the middle layer as $A$ and $B$ with site symmetry group $D_{3h}$, and the remaining positions as $A^{\prime}\mu$ and $B^{\prime}\mu$ with $\mu=1$ ($\mu=2$) for the top (bottom) layer with site symmetry group $C_{3v}$. The points $A$ and $B$ as well as the center of the hexagon in the middle layer are Wyckoff positions with multiplicity $m=1$, which we use as the origin of the three coordinate systems defined in Fig. 7. We therefore associate with each Wyckoff letter $W$ an index $(\alpha)$ corresponding to the coordinate systems $\alpha=a,b,c$. Ignoring the $z$ component, the positions of the C atoms in the unit cell for the three coordinate systems are

[TABLE]

where the superscripts $\alpha=a,b,c$ denote the coordinate system. The positions of the six atoms in unit cell $j$ are denoted by ${\bm{\mathrm{R}}}_{j}^{A}$, ${\bm{\mathrm{R}}}_{j}^{B}$, ${\bm{\mathrm{R}}}_{j}^{A^{\prime}\mu}$ and ${\bm{\mathrm{R}}}_{j}^{B^{\prime}\mu}$. In standard notation [58], the positions $A$, $B$ , $A^{\prime}\mu$ and $B^{\prime}\mu$ are characterized by the Wyckoff letters $c$, $e$, $g$, and $h$ respectively in coordinate system $(a)$, by $e$, $a$, $h$, and $i$ in coordinate system $(b)$, and by $a$, $c$, $i$, and $g$ in coordinate system $(c)$, see Table 1.

The point group $D_{3h}$ of the crystal structure of TLG is the same as for MoS2 with coordinate system shown in Fig. 2(a). The group of the wave vector at ${\bm{\mathrm{K}}}$ is $\mathcal{G}_{{\bm{\mathrm{K}}}}=C_{3h}$ with coordinate system shown in Fig. 2(b). At the three points ${\bm{\mathrm{M}}}_{i}$ the group of the wave vector is $\mathcal{G}_{{\bm{\mathrm{M}}}}=C_{2v}$ with coordinate system shown in Figs. 2(c)–2(e), which contains the twofold rotation $C_{2}$ about the dashed line, reflection $\sigma_{v}^{\prime}$ about the dashed line and the $z$ axis, and reflection $\sigma_{v}$ about the $xy$ plane.

IV.2 IRs of Bloch states in graphene

In graphene, the bands near the Fermi level are dominated by the $p$ orbitals of the C atoms. For SLG, using the coordinate systems in Fig. 4, the IRs of the symmetry-adapted $p$ orbitals are listed in Table 11 for the points ${\bm{\mathrm{\Gamma}}}$, ${\bm{\mathrm{K}}}$, and ${\bm{\mathrm{M}}}$ with groups of the wave vector $D_{6h}$, $D_{3h}$, and $D_{2h}$, respectively. For BLG, Fig. 6 shows the coordinate systems used for the points ${\bm{\mathrm{\Gamma}}}$, ${\bm{\mathrm{K}}}$, and ${\bm{\mathrm{M}}}$ with group of the wave vector $D_{3d}$, $D_{3}$, and $C_{2h}$, respectively. The IRs of the $p$ orbitals at these points are listed in Table 11. The coordinate systems used for TLG are the same as the ones for MoS2 (Fig. 2), so that the IRs of the $p$ orbitals can be taken from Table 4. The symmetry-adapted plane waves with the corresponding IRs are summarized in Table 6.

The IRs of the full Bloch functions for SLG, BLG, and TLG are listed in Table 12. The sets of Bloch states transforming as an IR in $\mathcal{G}_{{\bm{\mathrm{k}}}}$ for the wave vectors ${\bm{\mathrm{k}}}={\bm{\mathrm{\Gamma}}},{\bm{\mathrm{K}}},{\bm{\mathrm{M}}}$ are listed in Tables 13, 14, and 15 for SLG, Tables 16, 17, and 18 for BLG, and Tables 19, 20, and 21 for TLG.

V Selection Rules: Effect of Band IR Rearrangement

We show in the following that the selection rules for the observable matrix elements of a Hermitian operator $\mathcal{O}$ taken between Bloch states are not affected by the rearrangement of band IRs discussed in Sec. II.7 provided the perturbation of a crystal represented by the operator $\mathcal{O}$ preserves translational invariance. This condition for the operator $\mathcal{O}$ is certainly obeyed by the dipole operator $\mathcal{O}={\bm{\mathrm{r}}}$ representing optical transitions, but it does not apply to localized perturbations such as point defects for which anyway the specific location of the defect plays a crucial role [3]. We consider a translation ${\bm{\mathrm{\tau}}}$ of the coordinate system as introduced in Sec. II.7, so that $T=\{\openone\mskip 1.19998mu|\mskip 1.19998mu{\bm{\mathrm{\tau}}}\}$ is a unitary operator with $T^{\dagger}=T^{-1}$. This transforms both the states and operators. A state $|{\psi_{J\beta}}\rangle$ in the old coordinate system transforming according to the IR $\Gamma_{J}$ becomes $|{\psi_{J^{\prime}\beta^{\prime}}}\rangle=T|{\psi_{J\beta}}\rangle$ transforming according to $\Gamma_{J^{\prime}}=\Gamma_{\bm{\mathrm{\tau}}}\times\Gamma_{J}$, with $\Gamma_{\bm{\mathrm{\tau}}}=\{\mathcal{D}_{\bm{\mathrm{\tau}}}^{\bm{\mathrm{k}}}(g)=\exp(-i{\bm{\mathrm{b}}}_{g}\cdot{\bm{\mathrm{\tau}}}):g\in\mathcal{G}_{{\bm{\mathrm{k}}}}\}$ [Eq. (67)]. Similarly, $\langle{\psi_{I\alpha}}|$ transforming according to the IR $\Gamma_{I}^{\ast}$ becomes $\langle{\psi_{I^{\prime}\alpha^{\prime}}}|=\langle{\psi_{I\alpha}}|T^{\dagger}$ transforming according to $\Gamma_{I^{\prime}}^{\ast}=\Gamma_{\bm{\mathrm{\tau}}}^{\ast}\times\Gamma_{I}^{\ast}$. Invariance of the operator $\mathcal{O}$ under translations $T$ implies $[\mathcal{O},T]=0$, or

[TABLE]

so that both $\mathcal{O}$ and $\mathcal{O}^{\prime}$ transform according to the same representation (which need not be irreducible)

[TABLE]

Hence, we get (note $\Gamma_{\bm{\mathrm{\tau}}}^{\ast}\times\Gamma_{\bm{\mathrm{\tau}}}=\Gamma_{1}$)

[TABLE]

i.e., in the unprimed and primed coordinate system we get the same selection rules. In the primed coordinate system, the matrix elements become

[TABLE]

Note that for multidimensional IRs $\Gamma_{J}$ and $\Gamma_{J^{\prime}}=\Gamma_{\bm{\mathrm{\tau}}}\times\Gamma_{J}$, we can always choose the representation matrices such that $T|{\psi_{J\beta}}\rangle=|{\psi_{J^{\prime}\beta}}\rangle$, and similarly $\langle{\psi_{I\alpha}}|T^{\dagger}=\langle{\psi_{I^{\prime}\alpha}}|$, so that Eq. (105) becomes $\mathcal{O}_{\alpha\beta}^{\prime}=\mathcal{O}_{\alpha\beta}$.

VI Theory of Invariants

VI.1 General theory

The Bloch eigenstates $\Psi_{n{\bm{\mathrm{k}}}}^{I\beta}({\bm{\mathrm{r}}})$ at a wave vector ${\bm{\mathrm{k}}}$ transform according to an IR $\Gamma_{I}$ of the group of the wave vector $\mathcal{G}_{{\bm{\mathrm{k}}}}$. The knowledge of these IRs suffices to construct the general form of the effective Hamiltonian characterizing the Bloch eigenstates near the expansion point ${\bm{\mathrm{k}}}$ consistent with the symmetry operations in $\mathcal{G}_{{\bm{\mathrm{k}}}}$. This method is known as the theory of invariants [3]. The Hamiltonian can be expressed in terms of a general tensor operator denoted as ${\bm{\mathcal{K}}}$ that may depend on, e.g., the kinetic wave vector ${\bm{\mathrm{\kappa}}}$ measured from ${\bm{\mathrm{k}}}$, external electric and magnetic fields ${\bm{\mathcal{E}}}$ and ${\bm{\mathcal{B}}}$, strain ${\bm{\mathrm{\epsilon}}}$ and spin ${\bm{\mathrm{S}}}$. For all group elements $g\in\mathcal{G}_{{\bm{\mathrm{k}}}}$, the Hamiltonian $\mathcal{H}({\bm{\mathcal{K}}})$ obeys the invariance condition [3]

[TABLE]

According to the theory of invariants, each block $\mathcal{H}_{II^{\prime}}$ of the matrix $\mathcal{H}$ corresponding to a pair of bands transforming according to the IRs $\Gamma_{I}$ and $\Gamma_{I^{\prime}}$ has the form

[TABLE]

where $J$ labels the $L_{J}$-dimensional IRs $\Gamma_{J}$ contained in the product representation $\Gamma_{I}^{\ast}\times\Gamma_{I^{\prime}}$, the basis matrices $X_{l}^{J}$ and the irreducible tensor operators $\mathcal{K}_{J}^{i}$ constructed from the perturbations ${\bm{\mathcal{K}}}$ transform according to the IR $\Gamma_{J}$, and $a_{iJ}^{II^{\prime}}$ are constant prefactors. The index $i$ labels the irreducible tensor operators transforming as $\Gamma_{J}$. In general, we have multiple blocks $\mathcal{H}_{II^{\prime}}({\bm{\mathcal{K}}})$ corresponding to different bands $n$ and $n^{\prime}$ transforming according to $\Gamma_{I}$ and $\Gamma_{I^{\prime}}$. To simplify the notation, we drop these additional indices.

VI.2 Effect of band IR rearrangements

A translation of the coordinate system by ${\bm{\mathrm{\tau}}}$ changes the IR of an eigenfunction $\Psi_{n{\bm{\mathrm{k}}}}^{I}({\bm{\mathrm{r}}})$ from $\Gamma_{I}$ to $\Gamma_{J}=\Gamma_{\bm{\mathrm{\tau}}}\times\Gamma_{I}$, see Eq. (68). We have

[TABLE]

where we used $\Gamma_{\bm{\mathrm{\tau}}}^{\ast}\times\Gamma_{\bm{\mathrm{\tau}}}=\Gamma_{1}$ for one-dimensional IRs $\Gamma_{\bm{\mathrm{\tau}}}$. Therefore, the IRs contained in $\Gamma_{J}^{\ast}\times\Gamma_{J^{\prime}}$ are equal to the IRs contained in $\Gamma_{I}^{\ast}\times\Gamma_{I^{\prime}}$. This implies that translations ${\bm{\mathrm{\tau}}}$ of the coordinate system resulting in a rearrangement of the IRs assigned to the eigenfunctions do not affect the invariant expansion of the Hamiltonian $\mathcal{H}({\bm{\mathcal{K}}})$.

VI.3 Invariant Hamiltonian for MoS2

As an application of the theory of invariants, we consider monolayer MoS2. In this material, the lowest conduction and highest valence band are at the points ${\bm{\mathrm{K}}}$ and ${\bm{\mathrm{K}}}^{\prime}$ [29], hence we focus on these points (where $\mathcal{G}_{{\bm{\mathrm{K}}}}=C_{3h}$). Table 3 lists the mapping of axial (${\bm{\mathrm{A}}}$) and polar vectors (${\bm{\mathrm{P}}}$) under the relevant symmetry operations. This allows one to confirm the examples for basis functions listed for $C_{3h}$ in Table 27. Crystal momentum ${\bm{\mathrm{\kappa}}}$ and an electric field ${\bm{\mathcal{E}}}$ transform like polar vectors, whereas spin ${\bm{\mathrm{S}}}$ and a magnetic field ${\bm{\mathcal{B}}}$ transform like axial vectors. Hence we immediately obtain from Table 27 the lowest-order tensor operators listed in the second column of Table 22. In general, we obtain higher-order tensor operators using the Clebsch-Gordan coefficients that are tabulated in, e.g., Ref. [56]. However, this procedure is greatly simplified if all IRs of the relevant group are one-dimensional, which holds for $C_{3h}$. In such a case, the higher-order tensor operators can be constructed using the multiplication table for the IRs that is reproduced for $C_{3h}$ in Table 30. Irreducible tensor operators are generally not unique. Given two irreducible tensor operators $\mathcal{K}$ and $\mathcal{K}^{\prime}$ transforming according to the same IR $\Gamma_{I}$, any linear combinations of these tensors transforms likewise irreducibly according to $\Gamma_{I}$ [5]. We exploit this freedom to choose linear combinations of irreducible tensors such as $\kappa_{-}^{3}$ and $\kappa_{+}^{3}$ that have also a well-defined behavior under time reversal symmetry (see Sec. VI.4).

We can also consider the effects of strain. When stress deforms a crystalline solid, the symmetry of the system is altered which changes the energy spectrum of the material. Suppose under a deformation a point ${\bm{\mathrm{r}}}$ in a solid undergoes a displacement ${\bm{\mathrm{u}}}({\bm{\mathrm{r}}})$. For small homogeneous strain, the symmetric strain tensor is defined as ($i,j=x,y,z$) [59]

[TABLE]

However, considering MoS2 as a quasi-2D material, strain due to a perpendicular stress component is not relevant. The components $\epsilon_{ij}$ transform like the symmetrized products $\{\kappa_{i},\kappa_{j}\}$ [3] so that we get the lowest-order operators listed in the second column of Table 22 while mixed higher-order tensor operators are listed in the third column of Table 22.

Since the IRs of $C_{3h}$ are all one-dimensional, the Hamiltonian blocks (107) at the points ${\bm{\mathrm{K}}}$ and ${\bm{\mathrm{K}}}^{\prime}$ of MoS2 are one-dimensional. The $1\times 1$ basis matrices $X_{1}^{J}$ can be absorbed into the prefactors $a_{iJ}^{II^{\prime}}$. Hence, Eq. (107) can be simplified to

[TABLE]

The IRs corresponding to the eleven bands at ${\bm{\mathrm{K}}}$ and ${\bm{\mathrm{K}}}^{\prime}$ that are dominated by the Mo $d$ and S $p$ orbitals are listed in Table 23. Here we focus on constructing a generic $6\times 6$ Hamiltonian consisting of a sequence of bands transforming as $\Gamma_{1/3/2},\Gamma_{2/1/3},\Gamma_{3/2/1},\Gamma_{4/6/5},\Gamma_{5/4/6}$, and $\Gamma_{6/5/4}$, respectively (realized, e.g., by the bands $v_{1},c_{1},c_{3},c_{2},v_{3}$, and $v_{2}$; additional bands will replicate the behavior obtained for these bands). For these bands, the Hamiltonian matrix elements $\mathcal{H}^{{\bm{\mathrm{K}}}}({\bm{\mathcal{K}}})_{II^{\prime}}$ contain tensors $\mathcal{K}$ transforming as

[TABLE]

see Table 30 (Appendix 29).

VI.4 Time reversal

At the points ${\bm{\mathrm{K}}}$ and ${\bm{\mathrm{K}}}^{\prime}=-{\bm{\mathrm{K}}}$ with $\mathcal{G}_{{\bm{\mathrm{K}}}}=C_{3h}$, the Bloch functions $\Psi_{{\bm{\mathrm{K}}}}^{I}$ and the corresponding time reversed functions $\Theta\,\Psi_{{\bm{\mathrm{K}}}}^{I}$ are linearly dependent on each other [case $(a)$ according to Eq. (79), see also Table 9]. Also, the eigenfunctions at ${\bm{\mathrm{K}}}$ can be mapped onto the eigenfunctions at $-{\bm{\mathrm{K}}}$ by a vertical reflection $R=\sigma_{v}^{(2)}$ or a $180^{\circ}$ rotation $C_{2}^{\prime(2)}$, which is case (2) as defined in Eq. (80). [These two operations $R$ are elements of the group $D_{3h}$, the point group of MoS2, see Fig. 2(a).] Hence, for a Bloch state $\Psi_{{\bm{\mathrm{K}}}}^{I}$ of band $I$, the time reversed state $\Theta\,\Psi_{{\bm{\mathrm{K}}}}^{I}$ and the spatially transformed state $R\,\Psi_{{\bm{\mathrm{K}}}}^{I}$, with $R=\sigma_{v}^{(2)},C_{2}^{\prime(2)}$, obey the linear relation

[TABLE]

where $t_{R}^{I}$ is a phase factor (a unitary matrix if the dimensions of the IRs $\Gamma_{I}$ and $\Gamma_{I}^{\ast}$ was larger than one) that depends on the choice for the operation $R$. For the phase conventions used in Table 24, we have

[TABLE]

The matrix $\mathcal{H}({\bm{\mathcal{K}}})$ must then satisfy the additional condition

[TABLE]

where $\mathcal{t}_{R}$ is a diagonal matrix with elements $(\mathcal{t}_{R})_{II}=t_{R}^{I}$, $\zeta=+1$ ($-1$) for quantities that are even (odd) under time reversal such as ${\bm{\mathcal{E}}}$ and ${\bm{\mathrm{\epsilon}}}$ (${\bm{\mathrm{\kappa}}}$, ${\bm{\mathcal{B}}}$, and ${\bm{\mathrm{S}}}$), $\ast$ denotes complex conjugation and $t$ transposition. The components of polar vectors ${\bm{\mathrm{P}}}$ and axial vectors ${\bm{\mathrm{A}}}$ transform under $\sigma_{v}^{(2)}$ and $C_{2}^{\prime(2)}$ as follows (see Table 3)

[TABLE]

It follows from Eq. (114) whether off-diagonal terms in $\mathcal{H}({\bm{\mathcal{K}}})$ have real or imaginary prefactors. Tensor operators that give rise to invariants with real prefactors are marked in bold in Table 22. Furthermore, condition (114) provides a general criterion which terms are allowed by time-reversal symmetry on the diagonal of the Hamiltonian [when the tensor operators must transform according to the identity IR $\Gamma_{1}$, see Eq. (111)]. Here, our phase conventions imply that invariants (with real prefactors) can only be formed from tensor operators marked in bold in Table 22. Thus, for example, on the diagonal of the Hamiltonian the third-order trigonal term $\kappa_{+}^{3}+\kappa_{-}^{3}=2\kappa_{x}(\kappa_{x}^{2}-3\kappa_{y}^{2})$, as well as the field-dependent terms $i\kappa_{-}\mathcal{E}_{+}-i\kappa_{+}\mathcal{E}_{-}=2(\kappa_{y}\mathcal{E}_{x}-\kappa_{x}\mathcal{E}_{y})$ and $\mathcal{B}_{+}S_{-}+\mathcal{B}_{-}S_{+}=2(\mathcal{B}_{x}S_{x}+\mathcal{B}_{y}S_{y})$ are allowed by symmetry and thus present in the Hamiltonian. However, the terms $i\kappa_{-}^{3}-i\kappa_{+}^{3}=2\kappa_{y}(3\kappa_{x}^{2}-\kappa_{y}^{2})$, $\kappa_{+}\mathcal{E}_{-}+\kappa_{-}\mathcal{E}_{+}=2(\kappa_{x}\mathcal{E}_{x}+\kappa_{y}\mathcal{E}_{y})$ and $i\mathcal{B}_{-}S_{+}-i\mathcal{B}_{+}S_{-}=2(\mathcal{B}_{y}S_{x}-\mathcal{B}_{x}S_{y})$ are allowed by spatial symmetry; but these terms are forbidden by time-reversal symmetry and hence do not appear in the Hamiltonian. (They are allowed, though, as off-diagonal terms coupling different bands transforming according to the same IR, when these terms will have imaginary prefactors.)

Using the transformation $R$, one can derive the effective Hamiltonian for the valley ${\bm{\mathrm{K}}}^{\prime}$ using the transformation

[TABLE]

Alternatively, time reversal can be used, see Eq. (114). Note that the definition of $\mathcal{H}^{{\bm{\mathrm{K}}}^{\prime}}({\bm{\mathcal{K}}})$ relative to $\mathcal{H}^{{\bm{\mathrm{K}}}}({\bm{\mathcal{K}}})$ depends on the phase conventions used for the basis functions $\Psi_{{\bm{\mathrm{K}}}^{\prime}}^{I}$ relative to the phase conventions used for $\Psi_{{\bm{\mathrm{K}}}}^{I}$.

VI.5 Analysis of the invariant expansion

The invariant Hamiltonian at the point ${\bm{\mathrm{K}}}$ of the BZ becomes in lowest order of the wave vector ${\bm{\mathrm{\kappa}}}$ and spin-orbit coupling (spin ${\bm{\mathrm{S}}}$)

[TABLE]

with

[TABLE]

where $E_{I}$, $\gamma_{ij}$, and $\lambda_{ij}$ are material-dependent real parameters.

The lowest-order strain-dependent terms become

[TABLE]

with real parameters $\xi_{ij}$. The effect of electric and magnetic fields ${\bm{\mathcal{E}}}$ and ${\bm{\mathcal{B}}}$ can also be included in the Hamiltonian (117). For the electric field, we add on the diagonal a scalar potential $e{\bm{\mathcal{E}}}\cdot{\bm{\mathrm{r}}}$ and we replace crystal momentum by kinetic momentum $\hbar{\bm{\mathrm{\kappa}}}=-i\hbar{\bm{\mathrm{\nabla}}}+e{\bm{\mathrm{A}}}$, where ${\bm{\mathrm{A}}}$ is the vector potential due to the magnetic field. In addition, we may have terms that depend explicitly on the fields ${\bm{\mathcal{E}}}$ and ${\bm{\mathcal{B}}}$. The in-plane components $\mathcal{E}_{x}$ and $\mathcal{E}_{y}$ transform spatially like the wave vector components $\kappa_{x}$ and $\kappa_{y}$. However, ${\bm{\mathcal{E}}}$ is even under time reversal symmetry, whereas ${\bm{\mathrm{\kappa}}}$ is odd. In lowest order, the ${\bm{\mathcal{E}}}$-dependent terms become

[TABLE]

with imaginary prefactors $\eta_{ij}$. Since ${\bm{\mathcal{B}}}$ transforms in the same way as the spin ${\bm{\mathrm{S}}}$ both under spatial transformations and time reversal, the ${\bm{\mathcal{B}}}$-dependent terms are similar to Eq. (118u)

[TABLE]

with real parameters $\beta_{ij}$.

We can project the multiband Hamiltonian (117) on a subspace containing the bands of interest using quasidegenerate perturbation theory or Löwdin partitioning [5]. This gives rise to higher-order terms in the effective Hamiltonian that may include also mixed terms proportional to products of ${\bm{\mathrm{\kappa}}}$, ${\bm{\mathrm{S}}}$, ${\bm{\mathrm{\epsilon}}}$, ${\bm{\mathcal{E}}}$, and ${\bm{\mathcal{B}}}$. In particular, in the presence of a magnetic field ${\bm{\mathcal{B}}}$, the components of crystal momentum $\hbar{\bm{\mathrm{\kappa}}}=-i\hbar{\bm{\mathrm{\nabla}}}+e{\bm{\mathrm{A}}}$ do not commute, $\left[\kappa_{x},\kappa_{y}\right]=(e/i\hbar)\mathcal{B}_{z}$, so that antisymmetrized products of $\kappa_{x}$ and $\kappa_{y}$ appearing in higher-order perturbation theory give rise to terms proportional to $\mathcal{B}_{z}$. Similarly, terms proportional to in-plane electric fields $\mathcal{E}_{x},\mathcal{E}_{y}$ appear because of the presence of the scalar potential $e{\bm{\mathcal{E}}}\cdot{\bm{\mathrm{r}}}$ and the relation $\left[r_{i},\kappa_{j}\right]=i\delta_{ij}$.

VII Conclusions

Starting from a TB approach, we have developed a comprehensive theory to derive the IRs characterizing the Bloch eigenstates in a crystal by decomposing the TB basis functions into localized symmetry-adapted atomic orbitals and crystal-periodic symmetry-adapted plane waves. Both the symmetry-adapted atomic orbitals and the symmetry-adapted plane waves can easily be tabulated, thus accelerating the design and exploration of new materials. The symmetry-adapted basis functions block-diagonalize the TB Hamiltonian, which naturally facilitates further analysis of the band structure. While our analysis focused for clarity on symmorphic space groups, our theory can readily be generalized to nonsymmorphic groups.

The present work was motivated by the goal to develop a systematic theory of effective multiband Hamiltonians for the dynamics of Bloch electrons in external fields that break the symmetry of the crystal structure. Yet our general symmetry analysis of Bloch states will likely be useful for other applications, too.

Acknowledgements.

We appreciate stimulating discussions with G. Burkard, A. Kormányos, and U. Zülicke This work was supported by the NSF under Grant No. DMR-1310199. Work at Argonne was supported by DOE BES under Contract No. DE-AC02-06CH11357.

Appendix A Projection operators

A general method to identify the IRs of the eigenstates of a Hamiltonian uses projection operators [19, 38, 3, 39]. Given a symmetry group $G$ with IRs $\Gamma_{I}$, we can project a general function $f({\bm{\mathrm{r}}})$ onto its components $f_{I\beta}({\bm{\mathrm{r}}})$ transforming according to the $\beta$th component of the IR $\Gamma_{I}$ of $G$. Here the projection operators $\Pi_{I\beta}$ are

[TABLE]

where $h$ is the order of the group $G$, $n_{I}$ is the dimensionality of the IR $\Gamma_{I}$, $\mathcal{D}_{I}(g)$ are the representation matrices for $\Gamma_{I}$, and $P(g)$ are the symmetry operators corresponding to $g\in G$. Often, we denote $P(g)$ simply as $g$. If we are not interested in a particular component $\beta$ of $\Gamma_{I}$, we can use the “coarse-grained” projection operators

[TABLE]

where $\chi_{I}(g)\equiv\operatorname{tr}\mathcal{D}_{I}(g)$ are the characters for $\Gamma_{I}$. For one-dimensional IRs, the operators $\Pi_{I}$ become equivalent to $\Pi_{I\beta}$. The projection operators obey the completeness relation

[TABLE]

Appendix B Rearrangement lemma for IRs

Generally, if $\Gamma_{0}$ is a one-dimensional IR of a group $G$, then for any IR $\Gamma_{I}$ of $G$, the product representation $\Gamma_{I^{\prime}}=\Gamma_{0}\times\Gamma_{I}$ is irreducible. If two IRs $\Gamma_{I}$ and $\Gamma_{J}$ of $G$ are (in)equivalent, then the IRs $\Gamma_{I^{\prime}}=\Gamma_{0}\times\Gamma_{I}$ and $\Gamma_{J^{\prime}}=\Gamma_{0}\times\Gamma_{J}$ are also (in)equivalent. Hence, a multiplication of the IRs of $G$ by $\Gamma_{0}$ simply rearranges the sequence of IRs of $G$. These statements can be proven as follows: For the element $g\in G$, we denote the characters for the IRs $\Gamma_{0}$ and $\Gamma_{I}$ by $\chi_{0}(g)$ and $\chi_{I}(g)$, respectively. Hence the character of the representation $\Gamma_{I^{\prime}}$ becomes $\chi_{I^{\prime}}(g)=\chi_{0}(g)\,\chi_{I}(g)$. Since $\Gamma_{0}$ is one-dimensional, its characters are also its unitary representation matrices obeying $|\chi_{0}(g)|^{2}=1$. Using the orthogonality relations for characters we get

[TABLE]

where $h$ is the order of the group. Hence $\Gamma_{I^{\prime}}=\Gamma_{0}\times\Gamma_{I}$ is irreducible if $\Gamma_{I}$ is irreducible. Now consider two IRs $\Gamma_{I^{\prime}}=\Gamma_{0}\times\Gamma_{I}$ and $\Gamma_{J^{\prime}}=\Gamma_{0}\times\Gamma_{J}$. Using again the orthogonality relations for characters we get

[TABLE]

so that $\Gamma_{I^{\prime}}$ and $\Gamma_{J^{\prime}}$ are indeed (in)equivalent if $\Gamma_{I}$ and $\Gamma_{J}$ are (in)equivalent.

Appendix C Group tables

Character tables and basis functions for the groups $D_{3h}$, $D_{2h}$, $C_{3h}$, $C_{2v}$, and $C_{2h}$ are reproduced in Tables 25 – 29 following the designations by Koster et al. [39]. The second column in each table gives the designations for the IRs used by Dresselhaus et al. [39]. Examples of basis functions that transform irreducibly according to the different IRs are listed in the last column. These functions are expressed in terms of the cartesian components of polar (P) and axial (A) vectors using the coordinate systems defined for the respective groups in Figs. 2, 4, and 6. Finally, we reproduce in Table 30 the multiplication table for the IRs of the group $C_{3h}$ [56].