of 8
Supplementary Information: Snowflake phononic topological insulator at the
nanoscale
I. FINITE ELEMENTS SIMULATIONS
The results shown in Fig. 1-3 of the main text have
been obtained by solving the eigenvalue problem
div
h
C
:
grad
+(
grad
)
T
i
=
2
%!
2
,
(1)
using the COMSOL finite element solver.
denotes
the complex three-dimensional wave function related to
the mechanical displacement field
u
=
Re
·
e
i
!
t
of
the crystal,
C
is the elasticity tensor, and
%
the mass
density. Here
:
is a short-hand for the tensor product,
[
C
: grad
]
ij
=
C
ijkl
@
l
k
.Forconcreteness,weenvi-
sioned a phononic crystal slab of thickness
220 nm
as in
state-of-the-art optomechanical devices [1].
II. FULL BAND STRUCTURES FOR VARIOUS
CENTRAL SNOWFLAKE RADII
Fig. 1 of the main text just shows the frequency range
of the band structure which is relevant for the emerging
edge states at the domain interface for the parameters
r
=0
and
r
= 200
nm (Fig. 1d,e of the main text).
To give a larger overview Fig. 1 displays the same band
structures again but in an extended frequency range. In
addition to that, band structures for various other cen-
tral snowflake radii
r
c
=
r
+
r
are depicted.
We are
especially interested in the size of the band gap separat-
ing the two Dirac cones. As expected, the finite element
simulations show that the band gap separating the Dirac
cones closes when all snowflakes have the same radius,
r
c
=
r
= 1800
nm. For larger central snowflake radii,
r
c
>r
(
r>
0
), the band gap increases monotonically
with the radius. In this case, a limit to the band gap
size is set by the geometry [
r
c
2500
nm because for
r
c
= 2500
nm the size of the bridge connecting the trian-
gles vanishes]. On the other hand, for smaller snowflake
radii,
r
c
<r
(
r<
0
), the band gap depends non-
monotonically on the radius and reaches a maximum for
r
200
nm (
r
c
= 1600
nm). We attribute this non-
monotonic behavior to the breakdown of the four-band
approximation.
III. ZIGZAG VS ARMCHAIR DOMAIN WALL
Edge states that are not topologically protected cru-
cially depend on the details of the boundary of the sys-
tem. For example, a graphene strip supports edge states
in the presence of zig-zag edges but not armchair edges.
In our present system, our e
ff
ective Hamiltonian has a
topological bulk and, thus, edge states will appear at
any smooth domain wall. In principle, one could nev-
ertheless expect a di
ff
erent behavior for sharp domain
walls. In Fig. 2, we show the band structure for two dif-
ferent strips supporting sharp domain walls of the zig-zag
type (panel a) and armchair type (panel b), respectively.
The domain walls are created by choosing opposite val-
ues of
r
in the two domains while keeping the six-fold
rotational symmetry within each unit cell. Comparing
the band structure for di
ff
erent domain wall configura-
tions, cf. (a) and (b), we observe a qualitatively di
ff
erent
behavior only for the bands depicted with faded color
which correspond to edge states localized on the physical
boundaries. This is not surprising because the mapping
to our e
ff
ective topological Hamiltonian Eq. (2) com-
pletely breaks down at a sharp interface with the vac-
uum. On the other hand, as far as the remaining bands
are concerned there is no clear di
ff
erence between the two
domain walls configurations. This includes the domain
wall bands (the bands with linear dispersion crossing in
the vicinity of the
-point). Most importantly, in both
cases no minigap could be resolved within the numerical
accuracy of our finite elements simulations. [Such a mini-
gap could be reasonably expected to appear due to the
sharp breaking of the six-fold symmetry at the domain
wall.]
IV. DERIVATION OF EFFECTIVE
HAMILTONIAN UP TO QUADRATIC TERMS
The e
ff
ective Dirac Hamiltonian Eq. (4) of the main
text is obtained by projecting the elasticity equations (1)
onto a fixed in-plane quasi-momentum
k
=(
k
x,
k
y
)
close
to the
point. The resulting Hamiltonian is further pro-
jected onto a Hilbert space spanned by the appropriate
set of four normal modes
|
,
,
k
i
labeled by the indices
=
±
1
,
=
±
1
. In position space, we have
,
,
k
(
x
,z
)=
e
i
k
·
r
,
(
x
,z
)
,
where
r
=(
x, y
)
and
z
are the inplane and out-of-plane
coordinates, respectively. The states
,
(
x
,z
)
at the
point are chosen to be eigenstates of the three-fold rota-
tions, and of the discrete translational symmetry
ˆ
T
a
as-
sociated to the original unit cell (in the original Wigner-
Seitz cell their quasi-momenta lie at the
K
and
K
0
points
for
=1
and
=
1
,respectively)fullfillingEq.(1)of
the main text. We note that the mirror symmetry is
not strictly speaking necessary. In a crystal where this
symmetry is not present, one could use the alternative
definition
|
,
i
=
ˆ
T
|
,
i
=
ˆ
R
|
,
i
.
This alternative definition of the basis
|
,
i
has the dis-
advantage of not fixing the relative phase of the states
2
Helvetica-Fontsize: 16px
LaTex-Fonstize: 16px
LaTex-Fonstize: 20px
Helvetica-Fontsize: 20px
LaTex-Fonstize: 24px
Helvetica-Fontsize: 24px
Default Font-Size:
Picture produced by:
Name:
Snowflake/TopoInsulator/FEM/
2016_11_15_UnitcellFolding/
Plot3.m
TI_BandStructuresLoop.pdf
L
N
N
L
N
N
L
N
N
L
N
N
L
N
N
L
N
N
L
N
N
L
N
N
L
N
N
0.25
0.5
0.75
1
1.25
1.5
1.75
1.44
1.46
1.48
1400nm
1500nm
1600nm
1700nm
1800nm
1900nm
2000nm
2100nm
2200nm
15.4MHz
10.8MHz
12.4MHz
24.7MHz
35.2MHz
43.0MHz
frequency in GHz
quasi momentum
Figure 1. Band structures for various central radii
r
c
of the snowflake phononic crystal slab, along a path passing the high
symmetry points of the Brillouin zone (red, Fig. 1e of the main text). Here, the modes symmetric to the
xy
-plane are depicted
with darker colors. The lower row shows the region around the Dirac cones of interest, with the gray shaded area indicating
the complete band gap arising due to the breaking of the discrete translational symmetry
T
a
.
Helvetica-Fontsize: 16px
LaTex-Fonstize: 16px
LaTex-Fonstize: 20px
Helvetica-Fontsize: 20px
LaTex-Fonstize: 24px
Helvetica-Fontsize: 24px
Default Font-Size:
Picture produced by:
Name:
TI_Domainwall_Comparison.pdf
Snowflake/TopoInsulator/FEM/2016_11_17_StripeStraight/
PlotStraightVsArmchairVsZigZag.m
0
frequency in GHz
0
p
3
a
-
p
3
a
-
3
a
3
a
3
a
-
3
a
1
.
50
1
.
48
1
.
46
1
.
44
quasi momentum
p
3
a
·
k
x
a
b
Figure 2. Band structure of strip with zigzag edges (a) ver-
sus strip with armchair edges (b). Apart from the edge
states localized at the physical boundary of the system (de-
picted with faded colors) both band structures are identi-
cal. The simulated strips have the following dimensions:
(a)
r
l/u
=
200
nm,n
l
=18
,n
u
=17
and (b)
r
l/u
=
200
nm,n
l
=20
,n
u
=19
. See Fig. 3 of the main text for
concrete definitions of
n
u/l
.
For all FEM calculations fixed
boundaries (
u
=
0
)wereusedattheupperandlowerendof
the silicon slab.
with opposite
(because the time-reversal symmetry
ˆ
T
is anti-unitary). After fixing appropriately this rel-
ative phase, one arrives nevertheless at the same e
ff
ec-
tive Hamiltonian Eq. (2), see Ref. [2] for the case
g
=0
.Adi
ff
erentchoiceoftherelativephasecorre-
sponds to a gauge transformation that will change the
e
ff
ective Hamiltonian Eq. (2) but, obviously, not the
physical properties, e. g. the presence of edge states or
the Kramers degeneracy.
The e
ff
ective Hamiltonian, including also the quadratic
terms, is most conveniently derived using a di
ff
erent set
of eigenstates,
|
p
+
k
i
=
1
p
2
(
|
1
,
1
,
k
i
+
|
1
,
1
,
k
i
)
,
|
p
k
i
=
1
p
2
(
|
1
,
1
,
k
i
|
1
,
1
,
k
i
)
,
|
d
+
k
i
=
1
p
2
(
|
1
,
1
,
k
i
+
|
1
,
1
,
k
i
)
,
|
d
k
i
=
1
p
2
(
|
1
,
1
,
k
i
+
|
1
,
1
,
k
i
)
.
(2)
At the
point, the states
|
p
±
k
=
0
i
and
|
d
±
k
=
0
i
are of the
p
-
and
d
-type as it can be readily verified by using
ˆ
R
/
3
=
3
ˆ
R
2
/
3
ˆ
R
where
ˆ
R
2
/
3
|
,
i
=exp[
i
2
/
3]
|
,
i
,
to
find the expected behavior for
p
-and
d
-type orbital
under
60
-degree rotations,
ˆ
R
/
3
|
p
±
i
=exp[
i
/
3]
|
p
±
i
,
ˆ
R
/
3
|
d
±
i
=exp[
2
i
/
3]
|
d
±
i
.Fromthedefinitionsof
the
ˆ
’s and
ˆ
’s matrices, it also readily follows that the
p
-and
d
-type orbitals are eigenstates of the gap-opening
operator
ˆ
x
and of the conserved pseudo-spin
ˆ
S
x
ˆ
z
,
|
d
±
i
=
|
d
±
i
,
ˆ
x
|
p
±
i
=
|
p
±
i
,
ˆ
S
|
p
±
i
=
|
p
±
i
,and
ˆ
S
|
d
±
i
=
|
d
±
i
.FromEq.(2)ofthemaintextone
can also deduce the behavior under reflections
ˆ
M
xz
|
p
±
i
=
|
p
i
,
ˆ
M
xz
|
d
±
i
=
|
d
i
.
Below, we denote the eigenstates
defined in (2) by the more compact alternative notation
|
m,
k
i
where
m
=
±
1
,
±
2
is the six-fold quasi-angular
momentum carried by the pseudo-spin of the phonons,
1
,
k
i
=
|
p
±
k
i
and
2
,
k
i
=
|
d
±
k
i
.
The e
ff
ective Hamiltonian can be expanded in powers
of the quasi-momentum,
ˆ
H
k
=
ˆ
H
(0)
k
+
ˆ
H
(1)
k
+
ˆ
H
(2)
k
+
O
(
k
3
)
,
(3)
where the matrix elements of
ˆ
H
(0)
k
are independent of
k
,
those of
ˆ
H
(1)
k
are linear in
k
,andthoseof
ˆ
H
(2)
k
quadratic
in
k
. Next, we will show how each term
ˆ
H
(
i
)
k
,
i
=0
,
1
,
2
,
is constrained by the symmetries of the problem. In the
proof below, we will also assume the additional condition
|
,
i
=
ˆ
T
|
,
i
, which, when plugged into Eq. (2),
yields
ˆ
T
|
m,
k
i
=
|
m,
k
i
for the eigenstates of the
60
-degree rotations. We note that this condition sets a
constraint on the sum of the phases of the two eigenstates
|
,
i
and
|
,
i
.Sinceallrelativephasesfortheba-
sis
|
,
i
were already fixed by Eq. (1) of the main text,
the additional condition corresponds to fix a global phase
for the basis
|
,
i
.
This will clearly not change the form
of the final Hamiltonian but makes the proof less cum-
bersome. Thus, the additional condition does not lead to
any loss of generality of the proof.
We start from
ˆ
H
(0)
k
whose most general form is
ˆ
H
(0)
k
=
X
mm
0
h
(0)
m,m
0
|
m,
k
ih
m
0
,
k
|
(4)
where
h
(0)
mm
0
=(
h
(0)
m
0
m
)
to ensure that
ˆ
H
(0)
k
is hermitian.
We enforce the rotational symmetry by requiring
X
k
ˆ
H
k
=
X
k
ˆ
R
/
3
ˆ
H
k
ˆ
R
/
3
,
(5)
which should hold at each order in
k
.
Taking into account
that
ˆ
R
/
3
|
m,
k
i
=exp[
im
/
3]
|
m,
ˆ
R
k
i
where
ˆ
R
k
is the
rotated quasi-momentum,
ˆ
R
=
1
2
1
p
3
p
31
,
we find
X
k
ˆ
R
/
3
ˆ
H
(0)
k
ˆ
R
/
3
=
X
k
X
mm
0
e
i
(
m
0
m
)
/
3
h
(0)
m,m
0
|
m,
ˆ
R
k
ih
m
0
,
ˆ
R
k
|
=
X
k
X
mm
0
e
i
(
m
0
m
)
/
3
h
(0)
m,m
0
|
m,
k
ih
m
0
,
k
|
=
X
k
X
mm
0
h
(0)
m,m
0
|
m,
k
ih
m
0
,
k
|
The last equality follows from the six-fold rotational sym-
metry Eq. (5) and it implies that
h
(0)
mm
0
=0
for
m
6
=
m
0
.
In other words, the states carrying di
ff
erent quasi-angular
momentum are not coupled in
ˆ
H
(0)
k
.Wecanfurthercon-
strain the matrix elements
h
(0)
mm
0
by taking advantage
of the time-reversal symmetry. Under the time-reversal
symmetry we have
X
k
ˆ
T
ˆ
H
(0)
k
ˆ
T
=
X
k
m
h
(0)
m,m
|
m,
k
ih
m,
k
|
=
X
k
m
h
(0)
m,m
|
m,
k
ih
m,
k
|
.
Thus, from
P
k
ˆ
T
ˆ
H
(0)
k
ˆ
T
=
P
k
ˆ
H
(0)
k
it follows that
ˆ
h
(0)
mm
=
ˆ
h
(0)
m
m
.
No further constraint follows from the
mirror symmetries. In order to arrive at Eq. (2) of the
main text, we also have to express the allowed terms
in terms of the
ˆ
’s and
ˆ
’s Pauli matrices. This is
readily done by rewriting Eq. (4) with the constraint
ˆ
h
(0)
mm
=
ˆ
h
(0)
m
m
in terms of the basis
|
,
,
k
i
using Eq.
(2). In this basis, we find
ˆ
H
(0)
k
=
!
+
g
ˆ
x
(6)
where
!
=(
ˆ
h
(0)
11
+
ˆ
h
(0)
22
)
/
2
is the frequency of the tip
of the cones for the unperturbed lattice
r
=0
,and
g
=(
ˆ
h
(0)
22
ˆ
h
(0)
11
)
/
2
is the mass.
Next we follow a similar procedure with
ˆ
H
(1)
k
whose
most general form is
ˆ
H
(1)
k
=
X
mm
0
k
·
h
(1)
m,m
0
|
m,
k
ih
m
0
,
k
|
.
Here, for each fixed
m
and
m
0
,
m, m
0
=
±
2
,
±
1
,
h
(1)
mm
0
is
a
2
-vector, and
h
(1)
mm
0
=(
h
(1)
m
0
m
)
to ensure that
ˆ
H
(1)
k
is
hermitian. Under a rotation by
60
degrees we find
X
k
ˆ
R
/
3
ˆ
H
(1)
k
ˆ
R
/
3
=
X
k
X
mm
0
k
·
h
(1)
m,m
0
e
i
(
m
0
m
)
/
3
|
m,
ˆ
R
k
ih
m
0
,
ˆ
R
k
|
=
X
k
X
mm
0
k
·
ˆ
R
h
(1)
m,m
0
e
i
(
m
0
m
)
/
3
|
m,
k
ih
m
0
,
k
|
=
X
k
X
mm
0
k
·
h
(1)
m,m
0
|
m,
k
ih
m
0
,
k
|
.
4
Thus, we find the constraint
ˆ
R
h
(1)
m,m
0
e
i
(
m
0
m
)
/
3
=
h
(1)
m,m
0
.
(7)
This implies that
h
(1)
mm
0
=0
,oritisaneigenvectorof
the matrix
ˆ
R
with eigenvalue
exp[
i
(
m
m
0
)
/
3]
.
Since
ˆ
R
has eigenvalues
exp[
±
i
/
3]
,
h
(1)
m,m
0
can be di
ff
erent from
zero only for
m
=
m
m
0
=
±
1
. In other words, the
linear terms in the quasi-momentum induce transitions
that change the quasi-angular momentum by one unit.
We note that the allowed transitions do not flip the spin
ˆ
S.
This underlies the block structure of the Dirac Hamil-
tonian (2) of the main text. By solving Eq. (7), we find
h
(1)
1
,
2
=(
h
(1)
2
,
1
)
=
h
(1)
+
1
i
,
h
(1)
2
,
1
=(
h
(1)
1
,
2
)
=
h
(1)
1
i
,
(8)
where
h
(1)
±
are two complex numbers. The Hamiltonian
is further constrained by the remaining symmetries. By
applying the time-reversal symmetry we find
X
k
ˆ
T
ˆ
H
(1)
k
T
=
X
k
X
mm
0
k
·
(
h
(1)
m,m
0
)
|
m,
k
ih
m
0
,
k
|
=
X
k
X
mm
0
k
·
(
h
(1)
m,
m
0
)
|
m,
k
ih
m
0
,
k
|
=
X
k
X
mm
0
k
·
h
(1)
m,m
0
|
m,
k
ih
m
0
,
k
|
.
(9)
From the last equality (that enforces the
ˆ
T
symmetry)
h
(1)
m,m
0
=(
h
(1)
m,
m
0
)
=
h
(1)
m
0
,
m
.
Using Eqs. (8), it follows
h
(1)
+
=
h
(1)
.
(10)
We can constrain even more
h
(1)
mm
0
by using the mirror
symmetry
ˆ
M
xz
.Wefind,
X
k
ˆ
M
xz
ˆ
H
(1)
k
ˆ
M
xz
=
X
k
X
mm
0
(
)
m
+
m
0
k
·
h
(1)
m,m
0
|
m,
ˆ
M
k
ih
m
0
,
ˆ
M
k
|
=
X
k
X
mm
0
(
)
m
+
m
0
k
·
ˆ
M
h
(1)
m,
m
0
|
m,
k
ih
m
0
,
k
|
=
X
k
X
mm
0
k
·
h
(1)
m,m
0
|
m,
k
ih
m
0
,
k
|
,
(11)
where
ˆ
M
=
10
0
1
.
From the last equality in Eq. (11), we find
h
(1)
1
,
2
=
ˆ
M
h
(1)
1
,
2
,
and from Eq. (8)
h
(1)
+
=
h
(1)
.
From the above equation and Eq. (10), it follows that
h
(1)
+
=
h
(1)
and both matrix elements are real. As al-
ready outlined for the term
ˆ
H
(0)
k
,weexpressalso
ˆ
H
(1)
k
in
terms of the
’s and
’s Pauli matrices to obtain
ˆ
H
(1)
k
=
v
z
(
k
x
ˆ
x
k
y
ˆ
y
)
,
(12)
where
v
=
h
(1)
+
.
Thus, putting together Eqs. (3,6,12)
we find Eq. (2) of the main text.
Next, we go a step further, calculating the term
H
(2)
containing the quadratic contributions,
H
(2)
k
=
X
mm
0
(
k
ˆ
h
(2)
m,m
0
k
)
|
m,
k
ih
m
0
,
k
|
where for each fixed
m
and
m
0
,
m, m
0
2
±
2
,
±
1
,
ˆ
h
(2)
m,m
0
is
a
2
2
complex matrix that can be chosen to be symmet-
ric. Moreover,
ˆ
h
(2)
m,m
0
=(
ˆ
h
(2)
m
0
,m
)
to ensure that
H
(2)
k
is
hermitian. By applying a rotation by
60
-degrees we find
X
k
ˆ
R
(
/
3)
H
(2)
k
ˆ
R
(
/
3)
=
X
k
X
mm
0
k
ˆ
h
(2)
m,m
0
k
e
i
(
m
0
m
)
/
3
|
m,
ˆ
R
k
ih
m
0
,
ˆ
R
k
|
=
X
k
X
mm
0
(
ˆ
R
1
k
)
ˆ
h
(2)
m,m
0
(
ˆ
R
1
k
)
e
i
(
m
0
m
)
/
3
|
m,
k
ih
m
0
,
k
|
=
X
k
X
mm
0
k
ˆ
R
ˆ
h
(2)
m,m
0
ˆ
R
1
k
e
i
(
m
0
m
)
/
3
|
m,
k
ih
m
0
,
k
|
,
=
X
k
X
mm
0
k
ˆ
h
(2)
m,m
0
k
|
m,
k
ih
m
0
,
k
|
,
Using the last equality (that follows from the six-fold
rotational symmetry of
H
(2)
k
)
we find
e
i
(
m
0
m
)
/
3
ˆ
R
ˆ
h
(2)
m,m
0
ˆ
R
1
=
ˆ
h
(2)
m,m
0
.
(13)
For each fixed
m
and
m
0
,
m, m
0
2
±
2
,
±
1
,thisisjust
ahomogeneoussystemoffourlinearequationsoffour
variables (the entries of the matrix
ˆ
h
(2)
mm
0
). For the com-
binations of
m
and
m
0
where all equations are indepen-
dent, it has only a trivial solution
ˆ
h
(2)
m,m
0
=0
.Thus,
afinitematrixelement
ˆ
h
(2)
m,m
0
6
=0
is allowed only if
the system of linear equations defined in Eq. (13) is
not independent. It turns out that this is the case for
m
=(
m
0
m
)mod6 = 0
,
±
2
.Wenotethatthetran-
sitions that change the quasi-angular momentum by two
units flip the spin
ˆ
S
.Suchspin-flippingtransitionslift