Finding symmetries of the Kekule-Y continuum model¶
See also
The complete source code of this example can be found in kekule.py. A Jupyter notebook can be found in kekule.ipynb.
Kekule-Y¶
Find symmetries of Kekule-Y effective model, starting from the full hexagonal symmetry group with time reversal and particle-hole symmetry.
ham_kek_y = """
vs * kron((k_x * sigma_x + k_y * sigma_y), eye(2))
+ vt * kron(eye(2), (k_x * sigma_x + k_y * sigma_y))
"""
display(qsymm.sympify(ham_kek_y))
\[\begin{split}\displaystyle \left[\begin{matrix}0 & vt \left(k_{x} - i k_{y}\right) & vs \left(k_{x} - i k_{y}\right) & 0\\vt \left(k_{x} + i k_{y}\right) & 0 & 0 & vs \left(k_{x} - i k_{y}\right)\\vs \left(k_{x} + i k_{y}\right) & 0 & 0 & vt \left(k_{x} - i k_{y}\right)\\0 & vs \left(k_{x} + i k_{y}\right) & vt \left(k_{x} + i k_{y}\right) & 0\end{matrix}\right]\end{split}\]
H_kek_y= qsymm.Model(ham_kek_y, momenta=['k_x', 'k_y'])
candidates = qsymm.groups.hexagonal()
sgy, cg = qsymm.symmetries(H_kek_y, candidates)
Enumerate subgroups that protect double Dirac cone:
all_subgroups = qsymm.groups.generate_subgroups(sgy)
# Terms proportional to the identity matrix don't open a gap.
identity_terms = [qsymm.Model({qsymm.sympify(1): np.eye(4)})]
dd_sg = []
for subgroup, gens in all_subgroups.items():
family = qsymm.continuum_hamiltonian(gens, 2, 0)
family = qsymm.hamiltonian_generator.subtract_family(family, identity_terms)
if family == []:
dd_sg.append(subgroup)
dd_sg.sort(key=len)
dd_sg[0]
frozenset({1, R(π) C})
Twofold rotation combined with sublattice symmetry is enough. However 2-fold rotation is not a physical symmetry of the lattice model. Let’s get rid of 2-fold rotations:
rdd_sg = []
for subgroup, gens in all_subgroups.items():
invs = {
qsymm.inversion(2) * g
for g in [
qsymm.identity(2),
qsymm.time_reversal(2),
qsymm.particle_hole(2),
qsymm.chiral(2),
]
}
# Skip subgroups that have inversion
if subgroup & invs:
continue
family = qsymm.continuum_hamiltonian(gens, 2, 0)
family = qsymm.hamiltonian_generator.subtract_family(family, identity_terms)
if family == []:
rdd_sg.append(subgroup)
rdd_sg.sort(key=len)
rdd_sg[0]
frozenset({1, R(2π/3), R(-2π/3), C, R(2π/3) C, R(-2π/3) C})
Sublattice symmetry with 3-fold rotations is the minimal symmetry group required. Let’s test if any of the solutions has no sublattice symmetry:
any(g.antisymmetry for sg in rdd_sg for g in sg)
True
Sublattice symmetry is required for the protection of the double Dirac cone.
Continuous rotation symmetries of continuum models¶
Kekule-Y¶
_, cg = qsymm.symmetries(H_kek_y, candidates=[], continuous_rotations=True, prettify=True)
cg
\[\begin{split}e^{-i\phi L}H(\mathbf{k})e^{i\phi L} = H(R_{\phi}\mathbf{k}) \\R_{{\phi}} = R\left(\phi\right)\\L = \begin{bmatrix}1.+0.j&0.+0.j&0.+0.j&0.+0.j\\0.+0.j&0.+0.j&0.+0.j&0.+0.j\\0.+0.j&0.+0.j&0.+0.j&0.+0.j\\ 0.+0.j& 0.+0.j& 0.+0.j&-1.+0.j\\\end{bmatrix} \\\end{split}\]
Kekule-O¶
ham_kek_o= """
vs * (
k_x * kron(eye(2), sigma_x)
+ k_y * kron(sigma_z, sigma_y)
)
+ Delta * kron(sigma_y, sigma_y)
"""
display(qsymm.sympify(ham_kek_o))
\[\begin{split}\displaystyle \left[\begin{matrix}0 & vs \left(k_{x} - i k_{y}\right) & 0 & - \Delta\\vs \left(k_{x} + i k_{y}\right) & 0 & \Delta & 0\\0 & \Delta & 0 & vs \left(k_{x} + i k_{y}\right)\\- \Delta & 0 & vs \left(k_{x} - i k_{y}\right) & 0\end{matrix}\right]\end{split}\]
H_kek_o= qsymm.Model(ham_kek_o, momenta=['k_x', 'k_y'])
_, cg = qsymm.symmetries(H_kek_o, candidates=[], continuous_rotations=True, prettify=True)
cg
\[\begin{split}\left[ H(\mathbf{{k}}), L \right] = 0 \\L = \begin{bmatrix} 0. +0.j&-0. +0.j&-0. +0.j& 0.5+0.j\\-0. -0.j& 0. +0.j& 0.5-0.j&-0. -0.j\\-0. -0.j& 0.5+0.j& 0. +0.j&-0. +0.j\\ 0.5-0.j&-0. +0.j&-0. -0.j& 0. +0.j\\\end{bmatrix} \\\\e^{-i\phi L}H(\mathbf{k})e^{i\phi L} = H(R_{\phi}\mathbf{k}) \\R_{{\phi}} = R\left(\phi\right)\\L = \begin{bmatrix} 0.5+0.j&-0. -0.j&-0. -0.j&-0. +0.j\\-0. +0.j&-0.5-0.j& 0. -0.j&-0. +0.j\\-0. +0.j& 0. +0.j&-0.5+0.j&-0. +0.j\\-0. -0.j&-0. -0.j&-0. -0.j& 0.5+0.j\\\end{bmatrix} \\\end{split}\]
The rotation generators are different in the two cases
Lattice models¶
Generate Kekule-O model¶
6 sites (3 A and 3 B) per unit cell, spinless orbitals. Time-reversal, sublattice and 6-fold rotation and mirror symmetry.
sites = ['A1', 'B1', 'A2', 'B2', 'A3', 'B3']
norbs = [(site, 1) for site in sites]
# Time reversal
TR = qsymm.PointGroupElement(sympy.eye(2), True, False, np.eye(6))
# Chiral symmetry
C = qsymm.PointGroupElement(sympy.eye(2), False, True, np.kron(np.eye(3), ([[1, 0], [0, -1]])))
# Atom A rotates into B, B into A.
sphi = 2*sympy.pi/6
RC6 = sympy.Matrix([[sympy.cos(sphi), -sympy.sin(sphi)],
[sympy.sin(sphi), sympy.cos(sphi)]])
permC6 = {'A1': 'B1', 'B1': 'A2', 'A2': 'B2', 'B2': 'A3', 'A3': 'B3', 'B3': 'A1'}
C6 = qsymm.groups.symmetry_from_permutation(RC6, permC6, norbs)
RMx = sympy.Matrix([[-1, 0], [0, 1]])
permMx = {'A1': 'A3', 'B1': 'B2', 'A2': 'A2', 'B2': 'B1', 'A3': 'A1', 'B3': 'B3'}
Mx = qsymm.groups.symmetry_from_permutation(RMx, permMx, norbs)
symmetries = [C, TR, C6, Mx]
# Only need the unique hoppings, one weak and one strong
hopping_vectors = [
('A1', 'B1', np.array([0, -1])), # around ring
('A2', 'B3', np.array([0, -1])), # next UC
]
# construct Hamiltonian terms to symmetry-allowed terms
family = qsymm.bloch_family(hopping_vectors, symmetries, norbs=norbs)
qsymm.display_family(family)
\[\begin{split}\displaystyle \left[\begin{matrix}0 & e^{- i k_{y}} & 0 & 0 & 0 & e^{\frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}}\\e^{i k_{y}} & 0 & e^{- \frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}} & 0 & 0 & 0\\0 & e^{\frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}} & 0 & e^{\frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}} & 0 & 0\\0 & 0 & e^{- \frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}} & 0 & e^{i k_{y}} & 0\\0 & 0 & 0 & e^{- i k_{y}} & 0 & e^{\frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}}\\e^{- \frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}} & 0 & 0 & 0 & e^{- \frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}} & 0\end{matrix}\right]\end{split}\]
\[\begin{split}\displaystyle \left[\begin{matrix}0 & 0 & 0 & e^{\frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}} & 0 & 0\\0 & 0 & 0 & 0 & e^{- \frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}} & 0\\0 & 0 & 0 & 0 & 0 & e^{- i k_{y}}\\e^{- \frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}} & 0 & 0 & 0 & 0 & 0\\0 & e^{\frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}} & 0 & 0 & 0 & 0\\0 & 0 & e^{i k_{y}} & 0 & 0 & 0\end{matrix}\right]\end{split}\]
Make terms at k=0 and check degeneracy for some linear combination.
Gamma_terms = [term.subs(dict(k_x=0, k_y=0)) for term in family]
evals, evecs = np.linalg.eigh(
np.sum([
list(term.values())[0] * coeff
for term, coeff in zip(Gamma_terms, [1, 2])
],
axis=0)
)
evals
array([-4., -1., -1., 1., 1., 4.])
Find symmetry representation in the low energy subspace
# projector to low energy subspace
proj = evecs[:, 1:5]
UC3 = (C6 * C6).U
# projected symmetry operator
pC3 = proj.T.conj() @ UC3 @ proj
assert np.allclose(pC3 @ pC3.T.conj(), np.eye(4))
evals2, evecs2 = np.linalg.eig(pC3)
# Check angles of the eigenvalues
np.angle(evals2) / (2 * np.pi)
array([ 0.33333333, -0.33333333, 0.33333333, -0.33333333])
Generate Kekule-Y model¶
6 sites (3 A and 3 B) per unit cell, spinless orbitals. Time-reversal, sublattice and 3-fold rotation and mirror symmetry.
sites = ['A1', 'B1', 'A2', 'B2', 'A3', 'B3']
norbs = [(site, 1) for site in sites]
# Time reversal
TR = qsymm.time_reversal(2, U=np.eye(6))
# Chiral symmetry
C = qsymm.chiral(2, U=np.kron(np.eye(3), ([[1, 0], [0, -1]])))
# Atom A rotates into B, B into A.
sphi = 2*sympy.pi/3
RC3 = sympy.Matrix([
[sympy.cos(sphi), -sympy.sin(sphi)],
[sympy.sin(sphi), sympy.cos(sphi)],
])
permC3 = {'A1': 'A1', 'B1': 'B3', 'A2': 'A2', 'B2': 'B1', 'A3': 'A3', 'B3': 'B2'}
C3 = qsymm.groups.symmetry_from_permutation(RC3, permC3, norbs)
RMx = sympy.Matrix([[-1, 0], [0, 1]])
permMx = {'A1': 'A1', 'B1': 'B1', 'A2': 'A3', 'B2': 'B3', 'A3': 'A2', 'B3': 'B2'}
Mx = qsymm.groups.symmetry_from_permutation(RMx, permMx, norbs)
symmetries = [C, TR, C3, Mx]
# Only need the unique hoppings, one weak and one strong
hopping_vectors = [
('A1', 'B1', np.array([0, -1])), # Y
('A2', 'B3', np.array([0, -1])), # other
]
# construct Hamiltonian terms to symmetry-allowed terms
family = qsymm.bloch_family(hopping_vectors, symmetries, norbs=norbs)
qsymm.display_family(family)
\[\begin{split}\displaystyle \left[\begin{matrix}0 & e^{- i k_{y}} & 0 & e^{\frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}} & 0 & e^{\frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}}\\e^{i k_{y}} & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0\\e^{- \frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}} & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0\\e^{- \frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}} & 0 & 0 & 0 & 0 & 0\end{matrix}\right]\end{split}\]
\[\begin{split}\displaystyle \left[\begin{matrix}0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & e^{- \frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}} & 0 & e^{- \frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}} & 0\\0 & e^{\frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}} & 0 & e^{\frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}} & 0 & e^{- i k_{y}}\\0 & 0 & e^{- \frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}} & 0 & e^{i k_{y}} & 0\\0 & e^{\frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}} & 0 & e^{- i k_{y}} & 0 & e^{\frac{i k_{y}}{2}} e^{\frac{\sqrt{3} i k_{x}}{2}}\\0 & 0 & e^{i k_{y}} & 0 & e^{- \frac{i k_{y}}{2}} e^{- \frac{\sqrt{3} i k_{x}}{2}} & 0\end{matrix}\right]\end{split}\]
Make terms at k=0 and check degeneracy for some linear combination.
Gamma_terms = [term.subs(dict(k_x=0, k_y=0)) for term in family]
evals, evecs = np.linalg.eigh(
np.sum([
list(term.values())[0] * coeff
for term, coeff in zip(Gamma_terms, [1, 2])
],
axis=0)
)
evals
array([-5.19615242e+00, -8.69969719e-18, 0.00000000e+00, 0.00000000e+00,
0.00000000e+00, 5.19615242e+00])
Find symmetry representation of 3-fold rotations in the low energy subspace
# projector to low energy subspace
proj = evecs[:, 1:5]
UC3 = C3.U
# projected symmetry operator
pC3 = proj.T.conj() @ UC3 @ proj
assert np.allclose(pC3 @ pC3.T.conj(), np.eye(4))
evals2, evecs2 = np.linalg.eig(pC3)
# Check angles of the eigenvalues
np.angle(evals2) / (2 * np.pi)
array([ 0. , 0. , 0.33333333, -0.33333333])
The symmetry representation in the low energy subspace near the Gamma point is different in the two cases.