# Helical boundary modes from synthetic spin in a plasmonic lattice

Sang Hyun Park  Department of Electrical & Computer Engineering, University of
Minnesota, Minneapolis, Minnesota, 55455, USA  Michael Sammon  Department of
Electrical & Computer Engineering, University of Minnesota, Minneapolis,
Minnesota, 55455, USA  Eugene Mele  Department of Physics and Astronomy,
University of Pennsylvania, Philadelphia, Pennsylvania, 19104, USA  Tony Low
<EMAIL_ADDRESS>Department of Electrical & Computer Engineering, University of
Minnesota, Minneapolis, Minnesota, 55455, USA

###### Abstract

Artificial lattices have been used as a platform to extend the application of
topological physics beyond electronic systems. Here, using the two-dimensional
Lieb lattice as a prototypical example, we show that an array of disks which
each support localized plasmon modes give rise to an analog of the quantum
spin Hall state enforced by a synthetic time reversal symmetry. We find that
an effective next-nearest-neighbor coupling mechanism intrinsic to the
plasmonic disk array introduces a nontrivial $Z_{2}$ topological order and
gaps out the Bloch spectrum. A faithful mapping of the plasmonic system onto a
tight-binding model is developed and shown to capture its essential
topological signatures. Full wave numerical simulations of graphene disks
arranged in a Lieb lattice confirm the existence of propagating helical
boundary modes in the nontrivial band gap.

##  I Introduction

Artificial lattices can be patterned into ordered structures designed to
control efficient flow of energy and informationJoannapoulos2011 . There has
been particular interest in developing periodic structures for this purpose
that use nontrivial topology in their interior Bloch bandsLu2014a ;
Khanikaev2017 . These are known to support protected transport channels at
boundaries between distinct topological states in analogy to the well-studied
electronic surface states that occur at the boundaries of two and three
dimensional topological insulatorsHasan2010 . Indeed, nontrivial bulk Chern
bands in two dimensions have been successfully developed on photonicWang2009a
, plasmonicJin2017c  and even mechanical platformsNash2015 , where time
reversal symmetry is broken in order to establish their topological nature. In
striking contrast, photonic 2D analogs to quantum spin Hall states which
retain time reversal symmetry have been elusive. A key challenge is that the
Kramers degeneracies for half-integral angular momentum, whose connectivity is
essential for defining the topological state, are generically absent from
photonic and plasmonic analogsLu2016 . While multiple realizations of a
photonic quantum spin Hall state have been demonstratedHafezi2011 ;
Khanikaev2012 ; Wu2015 , intricate engineering via fine tuning of the design
is required to enforce the required symmetriesKhanikaev2017 .

In this work we suggest a robust resolution to this dilemma, and illustrate it
with studies on a two dimensional plasmonic Lieb latticeWeeks2010 ;
Goldman2011  as a prototype. Our model is characterized by two length (energy)
scales that express the confinement of excitations within and transmission
between nodes of a Lieb lattice network. Over a realizable parameter range in
which these scales are of comparable size, we identify both an emergent two-
fold Kramers degeneracies and a synthetic spin-orbit process that gaps the
Bloch spectra to endow the bands with nontrivial $Z_{2}$ topological
orderKane2005a . We test this idea with full wave numerical calculations to
confirm the existence of propagating helical modes on the boundaries and we
develop a mapping onto a tight binding model that captures the essential
spectral signature of this state. We emphasize that these phenomena in
artificial lattices are controlled by two energy scales which are tunable over
a wide dynamic range and access ordered states that are practically
unachievable in ordinary electronic materials. This strategy can be
generalized to a wide family of appropriately designed artificial lattices.

##  II Results and discussion

A plasmonic band structure can be calculated when the graphene nanodisks are
arranged in a periodic configuration. Here we consider placing graphene
nanodisks on a Lieb latticeGoldman2011 . When separation between the disks,
$a$, is much larger than the nanodisk diameter, $d$, coupling between the
plasmon modes is negligible and we observe a ladder of on-site resonant
multipole electromagnetic modes (see fig. 1a). The three groups of energy
bands each correspond to the dipole, quadrupole, and hexapole modes of the
graphene nanodisk. As the disk separation is decreased, coupling strength
between plasmon modes increases and the dispersion in the plasmon band
structure becomes stronger. In fig. 1b the plasmon band structure for
$a/d=1.05$ is shown. Interestingly, we find that the dipole, quadrupole, and
hexapole modes develop a non-trivial band structure resembling the band
structure of a tight binding Lieb lattice. Finally, when the disks begin to
overlap ($a<d$, fig. 1c), a stark change in the band structure is observed. A
low energy hydrodynamic plasmon mode for which the dispersion satisfies
$\omega\rightarrow 0$ as $q\rightarrow 0$ appears and the mapping onto a
tight-binding Lieb lattice breaks down. This situation corresponds to the
setup given in previous works Silveiro2013 ; Xiong2019 ; Jin2017c ; Jung2018
where periodic modulation of structure or density are used to perturb the
hydrodynamic mode.

To better understand the plasmon band structure given in fig. 1b, we develop a
mapping of the plasmonic crystal onto a tight-binding model. The governing
equations of the plasmonic response can be written as an eigenvalue equation
$\hat{H}\psi=\omega\psi$ whereJin2017c

$\hat{H}=\begin{pmatrix}0&\hat{V}\hat{\mathbf{p}}^{T}\sqrt{\omega_{F}(\mathbf{r})}\\\
\frac{e^{2}}{\pi\hbar}\sqrt{\omega_{F}(\mathbf{r})}\hat{\mathbf{p}}&0\end{pmatrix},\quad\psi=\begin{pmatrix}\Phi\\\
\mathbf{J}/\sqrt{\omega_{F}(\mathbf{r})}\end{pmatrix}.$ (1)

The eigenstate $\psi$ is a vector of the potential $\Phi$ and current density
$\mathbf{J}$, $\mathbf{\hat{p}}$ is the in-plane momentum operator, $\hat{V}$
is the Coulomb operator, and $\omega_{F}=E_{F}/\hbar$ is the Fermi frequency.
The localized plasmonic modes of the nanodisk are analogous to the atomic
orbitals in a tight-binding theory. From an eigenmode calculation of the
isolated graphene disk, we find that the electric potential is given by

$\Phi_{n,l}^{p}(\mathbf{r})=\begin{cases}u_{n}(r)\cos(l\theta),&p=+\\\
u_{n}(r)\sin(l\theta),&p=-\end{cases}$ (2)

where $u_{n}(r)$ gives the radial dependence of $\Phi$, $n$ is the number of
nodes in $u_{n}(r)$, and $l$ is a positive integer. Here we focus on the six
lowest energy modes which have $n=0$ and $l=1,2,3$. Henceforth the superscript
$n$ will be omitted for notational simplicity. Note that rotational symmetry
of the isolated graphene disk implies that $\psi_{l}^{+}$ and $\psi_{l}^{-}$
are degenerate in energy. Details on the Hamiltonian and eigenstates are given
in the supplementary information.

Equipped with the Hamiltonian and atomic orbitals, we may now calculate the
matrix elements of the tight-binding Hamiltonian. In general, the matrix
elements between orbitals $\psi_{l}^{p}$ and $\psi_{l^{\prime}}^{p^{\prime}}$
on disks separated by $\mathbf{R}$ are given by

$t_{ll^{\prime}}^{pp^{\prime}}(\mathbf{R})=-i\frac{e^{2}}{\pi\hbar}\int\Bigl{(}\mathbf{J}^{p*}_{l}(\mathbf{r})\cdot\nabla\Phi_{l^{\prime}}^{p^{\prime}}(\mathbf{r-R})+\Phi^{p*}_{l}(\mathbf{r})\nabla\cdot\mathbf{J}_{l^{\prime}}^{p^{\prime}}(\mathbf{r-R})\Bigr{)}d^{2}r$
(3)

where the inner product is defined such that $\hat{H}\psi=\omega\psi$ is a
Hermitian eigenproblemJin2017c  (see supplementary information for details).
Given that we have the current density and electric potential of the orbital
modes, it is straightforward to calculate the tight-binding model parameters.
For the Lieb lattice, the nearest neighbor (NN) hopping parameters can be
calculated by setting $\mathbf{R}=\pm
a\mathbf{\hat{x}}\equiv\pm\mathbf{a}_{x}$ or $\mathbf{R}=\pm
a\mathbf{\hat{y}}\equiv\pm\mathbf{a}_{y}$ in eq. 3.

We now take a closer look at the nearest neighbor coupling between the
quadrupole orbitals $\psi_{2}^{\pm}$. Either by explicit calculation through
eq. 3 or by inspection of the mode symmetries, it may be shown that hopping
from $\psi_{2}^{+}$ to $\psi_{2}^{-}$ given by $t_{22}^{+-}(\mathbf{a}_{j})$
is zero for $j=x,y$. It follows from the $C_{4}$ symmetry of the quadrupole
modes that hopping amplitudes in the $x$ and $y$ directions will be equal,
i.e. $t_{22}^{pp}(\mathbf{a}_{x})=t_{22}^{pp}(\mathbf{a}_{y})$ for $p=\pm$.
Finally, it is important to note that the $\psi_{2}^{+}$ and $\psi_{2}^{-}$
orbitals have hopping amplitudes with opposite sign. The plasmonic band
structure of the quadrupole with nearest neighbor hopping will therefore have
two copies of the the 3-band Lieb lattice band structure, where one of the
copies will be inverted in energy with respect to the other. This is indeed
what we observe for the quadrupole modes as shown in fig. 1b. A similar
analysis can be applied to the other multipole modes. Here we will be focusing
on the quadrupole mode band structure since it shares the $C_{4}$ symmetry of
the lattice.

We now extend the model to include orbital coupling between $\psi_{3}^{+}$ and
$\psi_{2}^{\pm}$. Note that the same analysis can be applied to coupling
between $\psi_{2}^{\pm}$ and $\psi_{1}^{\pm}$. In fig. 2a all allowed hoppings
between the $\psi_{3}^{+}$ and $\psi_{2}^{\pm}$ modes are shown along with
their respective signs. Hopping parameters $t_{23}^{++}(\pm\mathbf{a}_{y})$,
$t_{23}^{-+}(\pm\mathbf{a}_{x})$ are zero and are not shown. Applying a
$C_{4}$ rotation centered on the $\psi_{+3}$ mode transforms the configuration
of fig. 2a into fig. 2b with an overall sign change on all sites. This
transformation directly gives the hopping parameters between $\psi_{3}^{-}$
and $\psi_{2}^{\pm}$. All the allowed transitions are illustrated on an energy
level diagram in fig. 2c.

Considering second order hopping processes mediated by the $\psi_{3}^{\pm}$
orbital states, we find two processes that affect the quadrupole mode band
structure. First, a second order hopping along two links in the $x$ or $y$
directions (e.g. $t_{23}^{++}(\mathbf{a}_{x})$ followed by
$t_{32}^{++}(\mathbf{a}_{x})$) provides a uniform shift in the energy of all
orbital modes. A more interesting effect is found when a left or right turn is
made between the first and second hops. Taking into account the signs given in
fig. 2, we arrive at an effective second neighbor hopping given by

$W=it_{s}\sum_{pp^{\prime}}\sum_{\langle\langle
ij\rangle\rangle}v_{ij}c^{\dagger}_{ip}\left[\sigma_{y}\right]_{pp^{\prime}}c_{jp^{\prime}}$
(4)

where $v_{ij}=(\mathbf{d}_{1}\times\mathbf{d}_{2})_{z}$, $\mathbf{d}_{1,2}$
are the nearest neighbor bonds connecting site $i$ to $j$, $c^{\dagger}_{ip}$
is a creation operator for orbital $\psi_{2}^{p}$ on lattice site $i$, and
$\sigma_{y}$ is a Pauli matrix acting on the quadrupole orbital degrees of
freedom. This intrinsic interaction describing a coupling between the
$\psi_{2}^{\pm}$ orbitals that is path-dependent is similar to the Kane-Mele
spin-orbit interaction discussed in the context of electronic systemsKane2005a
. Strength of the interaction $t_{s}$ can be found by projecting the full
Hamiltonian into an effective Hamiltonian in the $\psi_{2}^{\pm}$ basis. The
effective Bloch Hamiltonian for the $\psi_{2}^{\pm}$ orbital basis may then be
written as

$h(\mathbf{k})=\begin{pmatrix}h_{+}(\mathbf{k})&w(\mathbf{k})\\\
w^{\dagger}(\mathbf{k})&h_{-}(\mathbf{k})\end{pmatrix}=\begin{pmatrix}h_{+}(\mathbf{k})&w(\mathbf{k})\\\
w^{\dagger}(\mathbf{k})&\alpha h_{+}(\mathbf{k})\end{pmatrix}$ (5)

where $h_{\pm}$ is the intra-orbital Hamiltonian for orbital $\psi_{2}^{\pm}$,
$\alpha$ is the ratio between hopping amplitudes $t_{22}^{++}$ and
$t_{22}^{--}$, and $w$ is a skew-Hermitian Hamiltonian that describes the
effective second neighbor coupling given by eq. 4. Matrix representations of
eq. 5 are given in the supplementary information. With the effective Bloch
Hamiltonian, we are able to map the quadrupole degrees of freedom of a
plasmonic lattice onto a simple tight-binding model with two basis orbitals on
each lattice site (see fig. 3a,b).

For the ideal case with $\alpha=-1$ in eq. 5, the Hamiltonian may be written
as

$h(\mathbf{k})=h_{+}(\mathbf{k})\otimes\sigma_{z}+w(\mathbf{k})\otimes
i\sigma_{y}.$ (6)

The band structure for this Hamiltonian is shown in fig. 3c. A strict
topological classification of this Hamiltonian based on its symmetries is
possible. Define a synthetic time-reversal operator as $T=\lambda_{z}\otimes
i\sigma_{y}K$ where $\lambda_{z}=\textrm{diag}(1,-1,1)$ is a sublattice
symmetry operator acting on the sublattice degrees of freedom, $\sigma_{y}$ is
a Pauli matrix acting on the orbital degrees of freedom, and $K$ represents
complex conjugation. Then, the Hamiltonian in eq. 6 satisfies
$Th(\mathbf{k})T^{-1}=h(-\mathbf{k})$. Since $T^{2}=-1$, a two-fold degeneracy
is enforced at the time reversal invariant momenta by Kramer’s theorem. A
particle-hole symmetry operator can also be defined as
$C=\lambda_{z}\otimes\sigma_{z}K$ such that
$Ch(\mathbf{k})C^{-1}=-h(-\mathbf{k})$ and $C^{2}=+1$. From the time-reversal
and particle-hole symmetry, we can now place the $\alpha=-1$ Hamiltonian into
the DIII symmetry classChiu2016 . In two-dimensions, this class is
characterized by a $Z_{2}$ topological invariant. Since inversion symmetry is
preserved, we may calculate the $Z_{2}$ invariantFu2007a  $\nu$ of a single
band by

$(-1)^{\nu}=\prod_{i}\delta_{i},\quad\delta_{i}=\xi(\Gamma_{i})$ (7)

where $\Gamma_{i}$ are the time reversal invariant momenta in the Brillouin
zone and $\xi(\Gamma_{i})$ is the parity eigenvalue of the band at
$\Gamma_{i}$. Unlike the electronic case for which $\delta_{i}$ is given as a
product of the parity eigenvalues for all occupied bands, here we associate a
separate $\nu$ with each band individually since there is no notion of band
filling. The $Z_{2}$ invariant for each band shown in fig. 3c reveals the
topologically non-trivial nature of the ideal plasmonic Lieb lattice.

Physically, $\alpha=-1$ requires fine tuning to a special state where the
symmetry is strictly enforced. The plasmonic lattice comprised of disks does
not have symmetries that enforce the $\alpha=-1$ condition. When
$\alpha\neq-1$, classification of the system into the DIII class no longer
holds and we do not have a sharp definition for the $Z_{2}$ topological
invariant. In addition, topologically trivial intra-orbital second neighbor
coupling is also present in the plasmonic lattice. Tight-binding calculations
for the bulk band structure with $\alpha=-1.2$ and the trivial second neighbor
coupling are shown in fig. 3d. Breaking of the synthetic time-reversal
symmetry is evident from the lifted degeneracy of the bulk bands at TRIM
points $\Gamma$ and $X$. The spectrum is also no longer symmetric with respect
to $E=0$ because particle-hole symmetry is broken by the trivial second
neighbor hopping.

However, by examining the edge states of the physical $\alpha\neq-1$, we find
that the topological signatures of the ideal lattice are still present. The
spectrum of the edge modes can be studied in a finite system ribbon geometry
and shows that the edge modes from the ideal DIII structure (fig. 4a) survive
even when the strict symmetry classification breaks down (fig. 4b). A
difference between the two cases is observed at the edge mode crossing points
shown in fig. 4d,e where the symmetry broken case shows an avoided crossing
originating from mixing of the edge modes. In principle, propagating edge
modes near this crossing point will be able to back-scatter. However, the size
of the edge mode avoided crossing is roughly two orders of magnitude smaller
than the band width, making it negligible in most practical situations. Full
wave numerical simulations of the bulk plasmonic band structure shown in fig.
3e confirm that while the $\alpha=-1$ condition is not strictly satisfied, the
edge modes of topological origin still survive as shown in fig. 4c,f.

To elucidate the properties of the edge modes, we apply a transformation into
the angular momentum basis to eq. 5. The angular momentum basis states are
given by $\tilde{\psi}_{m}=\psi_{|m|}^{+}+i\textrm{sgn}(m)\psi_{|m|}^{-}$
where the electric potential of $\tilde{\psi}_{m}$ is
$\tilde{\Phi}_{m}(\mathbf{r})\equiv u_{0}(r)e^{im\theta}$ and $m=\pm 1,\pm
2,\pm 3$. Applying the transformation to eq. 5 gives the tight binding
Hamiltonian in the angular momentum basis

$\tilde{h}(\mathbf{k})=\begin{pmatrix}\bar{h}(\mathbf{k})+iw(\mathbf{k})&\Delta(\mathbf{k})\\\
\Delta(\mathbf{k})&\bar{h}(\mathbf{k})-iw(\mathbf{k})\end{pmatrix}$ (8)

where $\bar{h}=(h_{+}+h_{-})/2$ and $\Delta=(h_{+}-h_{-})/2$. Written in this
form, it is clear that the plasmonic crystal will behave in a similar manner
to the electronic quantum spin-hall insulator where the electronic spin
degrees of freedom is mapped onto the chirality of the plasmon modes. The
plasmonic edge modes observed in fig. 4f will have a potential proportional to
$e^{\pm i2\theta}$ on each disk where the chirality is tied to the propagation
direction.

The existence of helical edge states can be confirmed by a simulation on a
finite array of graphene disks. In order to couple with the
$\tilde{\psi}_{-2}$ angular momentum mode, we use a circularly polarized
dipole with dipole moment $\mathbf{p}=\hat{\mathbf{x}}+i\hat{\mathbf{y}}$
placed at the edge of an A site graphene disk. Simulation results in fig. 5
show a projection of the excited plasmon fields onto the $\tilde{\psi}_{\pm
2}$ modes. We find that the dipole strongly couples to a $\tilde{\psi}_{-2}$
mode which propagates counter-clockwise along the boundaries of the system. A
weak excitation of $\tilde{\psi}_{+2}$ is caused by imperfect coupling between
the circularly polarized dipole and the $\tilde{\psi}_{-2}$ mode.

Finally, we identify two other perturbations that open a topologically trivial
gap at the M point and may occur in an array of plasmonic disks. The first is
an inversion breaking perturbation which may result from an asymmetric hopping
along either the $x$ or $y$ links in the lattice. A shift of the A or C site
plasmonic disk towards the B site will give rise to the inversion breaking
perturbation. Shifting the B site energy also opens up a gap on one side of
the flat bands. A shift in the Fermi energy of the plasmonic disk will result
in this type of perturbation. When applying either of the above mentioned
perturbations to an ideal DIII class system, the edge modes are removed only
after the gap closes and reopens. Interestingly, even in the physical case
where we are no longer in the ideal DIII class, the edge modes persist until
the gap is closed, implying that the edge modes remain robust to external
perturbations. Results on perturbations are given in the supplementary
information.

##  III Conclusion/Discussion

In summary, we have shown that propagating helical edge modes induced by an
intrinsic synthetic spin-orbit process exists in an artificial plasmonic
crystal. We consider the Lieb lattice as a prototypical example and develop a
mapping of the plasmonic system onto a simple tight-binding model. With the
tight-binding model, we identify an ideal limit of the system in which it can
be classified into the DIII symmetry class. In this limit, a non-trivial
$Z_{2}$ invariant is calculated and the associated topologically protected
edge states are shown to exist. It is important to note that this ideal limit
is distinct from the case of electron spins on a Lieb lattice and is a result
of the opposite sign hopping between the localized plasmon modes. Although the
physical plasmonic lattice does not strictly follow the behavior of this ideal
limit, we are able to verify that the system is close enough to retain the
edge states found for the ideal DIII limit. Propagation of the helical
plasmonic edge modes is verified using full wave numerical simulations.

Experimentally, the proposed helical edge modes can be excited by placing a
near-field antenna at the edge of the plasmon disk array as shown in fig. 5c.
For the setup used in this work, the gap in which the edge mode can be
observed is 4meV. We expect the edge mode in this gap to be observable as a
plasmonic band gap of similar size has been successfully resolved using the
near-field scanning microscope techniqueXiong2019 . The energy scale of the
plasmonic band structure can also be shifted by simply scaling the geometry of
the lattice (see supplementary information). For a disk size of $d=700$nm and
periodicity $a=735$nm, the quadrupole band energies are lowered to 70meV.
Importantly, the ratio between bandwidth and band gap stays constant as the
geometry is scaled. Hence the physical phenomena that we have examined for the
setup in this work can be translated to a different energy range by scaling
the geometry accordingly.

## Methods

All electromagnetic simulations were performed using COMSOL Multiphysics.
Graphene is modeled as a surface current with a Drude conductivity,
$\sigma(\omega)=\frac{e^{2}}{\pi\hbar}\frac{E_{F}}{\hbar\gamma-i\hbar\omega}$.
Scattering time is set to $\tau=1$ps where $\gamma=1/\tau$. The bulk plasmonic
band structure shown in fig. 1 and fig. 3 was performed by assuming periodic
boundary conditions along both the in-plane directions. In the out-of-plane
direction, the structure is padded with 1$\mu$m of free space and a perfectly
matched layer of thickness 300nm. The plasmonic band structure calculation for
a ribbon geometry shown in fig. 4c,f was performed using a ribbon length of 25
unit cells with periodic boundary conditions only in the direction parallel to
the ribbon edges. The circularly polarized dipole used for fig. 5c was placed
5nm above the edge of the graphene nanodisk. Projection of the potential onto
the angular momentum states $\tilde{\psi}_{\pm 2}$ was calculated as $\int
e^{i2\theta}\Phi(\mathbf{r})/|\Phi(\mathbf{r})|d^{2}r-\int
e^{-i2\theta}\Phi(\mathbf{r})/|\Phi(\mathbf{r})|d^{2}r$ for each disk.

###### Acknowledgements.

SHP, MS and TL acknowledge support by the National Science Foundation,
NSF/EFRI Grant No. EFRI-1741660. Work by EJM was supported by the Department
of Energy under grant DE-FG02-84ER45118. SHP acknowledges partial support by
the National Science Foundation through the University of Minnesota MRSEC
under Award Number DMR-2011401.

## References

  * (1) Joannapoulos, J., Johnson, S., Winn, J. & Meade, R.  _Photonic Crystals: Molding the Flow of Light_ (Princeton University Press, 2011), 2nd editio edn.
  * (2) Lu, L., Joannopoulos, J. D. & Soljačić, M.  Topological photonics.  _Nat. Photonics_ 8, 821–829 (2014).
  * (3) Khanikaev, A. B. & Shvets, G.  Two-dimensional topological photonics.  _Nat. Photonics_ 11, 763–773 (2017).
  * (4) Hasan, M. Z. & Kane, C. L.  Colloquium: Topological insulators.  _Rev. Mod. Phys._ 82, 3045–3067 (2010).
  * (5) Wang, Z., Chong, Y., Joannopoulos, J. D. & Soljačić, M.  Observation of unidirectional backscattering-immune topological electromagnetic states.  _Nature_ 461, 772–775 (2009).
  * (6) Jin, D. _et al._ Infrared Topological Plasmons in Graphene.  _Phys. Rev. Lett._ 118, 1–6 (2017).
  * (7) Nash, L. M. _et al._ Topological mechanics of gyroscopic metamaterials.  _Proc. Natl. Acad. Sci. U. S. A._ 112, 14495–14500 (2015).
  * (8) Lu, L., Joannopoulos, J. D. & Soljačić, M.  Topological states in photonic systems.  _Nat. Phys._ 12, 626–629 (2016).
  * (9) Hafezi, M., Demler, E. A., Lukin, M. D. & Taylor, J. M.  Robust optical delay lines with topological protection.  _Nat. Phys._ 7, 907–912 (2011).
  * (10) Khanikaev, A. B. _et al._ Photonic topological insulators.  _Nat. Mater._ 12, 233–239 (2012).
  * (11) Wu, L. H. & Hu, X.  Scheme for achieving a topological photonic crystal by using dielectric material.  _Phys. Rev. Lett._ 114, 1–5 (2015).
  * (12) Weeks, C. & Franz, M.  Topological insulators on the Lieb and perovskite lattices.  _Phys. Rev. B_ 82, 085310 (2010).
  * (13) Goldman, N., Urban, D. F. & Bercioux, D.  Topological phases for fermionic cold atoms on the Lieb lattice.  _Phys. Rev. A_ 83, 063601 (2011).
  * (14) Kane, C. L. & Mele, E. J.  Z2 Topological Order and the Quantum Spin Hall Effect.  _Phys. Rev. Lett._ 95, 146802 (2005).
  * (15) Silveiro, I., Manjavacas, A., Thongrattanasiri, S. & García de Abajo, F. J.  Plasmonic energy transfer in periodically doped graphene.  _New J. Phys._ 15, 033042 (2013).
  * (16) Xiong, L. _et al._ Photonic crystal for graphene plasmons.  _Nat. Commun._ 10, 4780 (2019).
  * (17) Jung, M., Fan, Z. & Shvets, G.  Midinfrared Plasmonic Valleytronics in Metagate-Tuned Graphene.  _Phys. Rev. Lett._ 121, 086807 (2018).
  * (18) Chiu, C.-K., Teo, J. C., Schnyder, A. P. & Ryu, S.  Classification of topological quantum matter with symmetries.  _Rev. Mod. Phys._ 88, 035005 (2016).
  * (19) Fu, L. & Kane, C. L.  Topological insulators with inversion symmetry.  _Phys. Rev. B_ 76, 045302 (2007).

Figure 1: Plasmonic band structure of graphene nanodisk arrays. a, Array
schematic and plasmonic band structure of graphene nanodisk array with disk
separation $a=200$nm and disk radius $d=100$nm. Bands corresponding to the
dipole(red), quadrupole(black), and hexapole(blue) are shown in the band
diagram. b, The disk separation is reduced to $a=105$nm while disk radius is
fixed to $d=100$nm. The dashed box indicates the unit cell. Bands formed from
the dipole, quadrupole, and hexapole modes of the graphene nanodisk with the
same color coding as a. c, For $a=85$nm and $d=100$nm the disks are now
overlapping. Color coding of the bands is not applied since classification in
terms of the multipole modes of the nanodisk does not hold in this case. For
all simulations, the graphene Fermi energy is set to $E_{F}=0.5$eV and the
substrate dielectric constant is $\epsilon_{s}=2.2$. Figure 2: Orbital
coupling between $\psi_{2}^{\pm}$ and $\psi_{3}^{\pm}$. a, Representation of
nearest neighbor hopping between $\psi_{3}^{+}$ and $\psi_{2}^{\pm}$. Hoppings
for $\mathbf{R}=\pm\mathbf{a}_{x},\pm\mathbf{a}_{y}$ are considered. The
electric potential $\Phi$ of the modes are shown. The plus and minus signs
associated with each hopping parameter indicate its sign. b, Hopping between
$\psi_{3}^{-}$ and $\psi_{2}^{\pm}$ are shown. All hopping parameters for b
can be generated by rotation of the configuration in a. c, The allowed hopping
transitions are represented on an energy level diagram. A second order
transition connecting $\psi_{2}^{+}$ to $\psi_{2}^{-}$ is possible via the
$\psi_{3}^{\pm}$ modes. Figure 3: Tight binding mapping of a plasmonic Lieb
lattice. a, Tight binding model onto which the plasmonic crystal can be mapped
onto. The red and blue arrows represent two paths through which the
interaction eq. 4 connects a $\psi_{2}^{+}$ state to a $\psi_{2}^{-}$ state.
b, Schematic representation of a graphene nanodisk Lieb lattice with
$\psi_{2}^{+}$ and $\psi_{2}^{-}$ orbitals on each site. Parameters are
equivalent to the setup in fig. 1b. c,d, Band structure calculated from tight
binding model given in eq. 5 with $\alpha=-1$ and $\alpha=-1.2$ respectively.
For the $\alpha=-1$ case, the $Z_{2}$ topological invariant $\nu$ for each
band is also shown. For the $\alpha=-1.2$ case, trivial second neighbor
hopping is included. e, Full wave electromagnetic simulation of plasmonic band
structure for configuration shown in b.  Figure 4: Edge modes of the plasmonic
Lieb lattice. a,d Tight-binding band structure calculation for Lieb lattice in
a ribbon geometry with $\alpha=-1$. The edge modes are highlighted in red
while the bulk band projections are colored in gray. All edge modes are doubly
degenerate. A more detailed view of the crossing point of a is given in d. b,e
Tight-binding band structure for $\alpha=-1.2$ and with trivial second
neighbor hopping included. c,f Full wave numerical simulation for plasmonic
Lieb lattice in ribbon geometry. Parameters are identical to those used for
fig. 1b and fig. 3e Figure 5: Dipole launching the edge states of the
plasmonic Lieb lattice.a, Band structure of the plasmonic Lieb lattice in a
ribbon geometry. Bands traversing the bulk gap are all doubly degenerate. b,
Schematic representation of the eigenstate corresponding to the black diamond
in a. The blue(red) shows projection onto state
$\tilde{\psi}_{-2}$($\tilde{\psi}_{+2}$). Arrows indicate the group velocity
of each edge mode. c, Simulation of topological edge plasmon launched by
circularly polarized dipole at frequency indicated in a. Magnitude of the
plasmon electric field is multiplied by a projection onto the
$\tilde{\psi}_{\pm 2}$ angular momentum states. Dipole source is placed at
black arrow 5nm above the graphene nanodisk.


# Utterance-level aggregation for speaker recognition in the wild

###### Abstract

The objective of this paper is speaker recognition ‘in the wild’ – where
utterances may be of variable length and also contain irrelevant signals.
Crucial elements in the design of deep networks for this task are the type of
trunk (frame level) network, and the method of temporal aggregation. We
propose a powerful speaker recognition deep network, using a ‘thin-ResNet’
trunk architecture, and a dictionary-based NetVLAD or GhostVLAD layer to
aggregate features across time, that can be trained end-to-end. We show that
our network achieves state of the art performance by a significant margin on
the VoxCeleb1 test set for speaker recognition, whilst requiring fewer
parameters than previous methods. We also investigate the effect of utterance
length on performance, and conclude that for ‘in the wild’ data, a longer
length is beneficial. ††http://www.robots.ox.ac.uk/ṽgg/research/speakerID

Index Terms—  speaker recognition, speaker verification, speech, deep
learning, CNNs

##  1 Introduction

Speaker recognition ‘in the wild’ has received an increasing amount of
interest recently due to the availability of free large-scale datasets [1, 2,
3], and the easy accessibility of deep learning frameworks [4, 5, 6]. For
speaker recognition, the goal is to condense information into a single
utterance-level representation, unlike speech recognition where frame-level
representations are desired. Obtaining a good utterance level representation
becomes particularly important for speech obtained under noisy and
unconstrained conditions, where irrelevant parts of the signal must be
filtered out. Therefore, a key area of research in deep learning for speaker
recognition is to investigate how to effectively aggregate frame-level
characteristics into utterance-level speaker representations.

Earlier deep neural network (DNN) based speaker recognition systems have
naïvely used pooling methods that have been successful for visual recognition
tasks, such as average pooling [2, 3, 7, 8] or fully connected layers [9, 10]
to condense frame-level information into utterance-level representations.
Although such methods serve the purpose of aggregating frame-level information
into a single representation whilst still allowing back-propagation, the
aggregation is not content dependent, so they are not able to consider which
parts of the input signal contain the most relevant information.

On the other hand, traditional methods for speaker and language identification
such as i-vector systems have explored the use of statistical or dictionary-
based methods for aggregation. A number of recent works have proposed to bring
similar methods to deep speaker recognition [11, 12, 13, 14, 15] (described in
Sec. 1.1). Based on these works, we propose to marry the best of both
Convolutional Neural Networks (henceforth, CNNs) and a dictionary-based
NetVLAD [16] layer, where the former is known for capturing local patterns,
and the latter can be discriminatively trained for aggregating information
into a fixed-sized descriptor from an input of arbitrary size, such that the
final representation of the utterance is unaffected by irrelevant information.

We make the following contributions: (i) We propose a powerful speaker
recognition deep network, based on a NetVLAD [16] or GhostVLAD [17] layer that
is used to aggregate ‘thin-ResNet’ architecture frame features; (ii) The
entire network is trained end-to-end using a large margin softmax loss on the
large-scale VoxCeleb2 [3] dataset, and achieves a significant improvement over
the current state-of-the-art verification performance on VoxCeleb1, despite
using fewer parameters than the current state-of-the-art architectures [3,
13]; and, (iii) We analyse the effect of input segment length on performance,
and conclude that for ‘in the wild’ sequences having longer utterances (4s or
more) is a significant improvement over shorter segments.

###  1.1 Related works

End-to-end deep learning based systems for speaker recognition usually follow
a similar three-stage pipeline: (i) frame level feature extraction using a
deep neural network (DNN); (ii) temporal aggregation of frame level features;
and (iii) optimisation of a classification loss. In the following, we review
the three components in turn.

The trunk DNN architecture used is often either a 2D CNN with convolutions in
both the time and frequency domain [2, 3, 15, 18, 19, 20], or a 1D CNN with
convolutions applied only to the time domain [11, 12, 13, 21]. A number of
papers [8, 22] have also used LSTM-based front-end architectures.

The output from the feature extractor is a variable length feature vector,
dependant on the length of the input utterance. Average pooling layers have
been used in [2, 3, 8] to aggregate frame-level feature vectors to obtain a
fixed length utterance-level embedding. [11] introduces an extension of the
method in which the standard deviation is used as well as the mean – this
method is termed statistical pooling, and used by [12, 21]. Unlike these
methods which ingest information from all frames with equal weighting, [20,
22] have employed attention models to assign weight to the more discriminative
frames. [13] combines the attention models and the statistical model to
propose attentive statistics pooling – this method holds the current state-of-
the-art performance on the VoxCeleb1 dataset. The final pooling strategy of
interest is the Learnable Dictionary Encoding (LDE) proposed by [14, 15]. This
method is closely based on the NetVLAD layer [16, 23] designed for image
retrieval.

Typically, such systems are trained end-to-end for classification with a
softmax loss [13] or one of its variants, such as the angular softmax [15]. In
some cases, the network is further trained for verification using the
contrastive loss [2, 3, 24] or other metric learning losses such as the
triplet loss [7]. Similarity metrics like the cosine similarity or PLDA are
often adopted to generate a final pairwise score.

##  2 Methods

For speaker recognition, the ideal model should have the following properties:
(1) It should ingest arbitrary time lengths as input, and produce a fixed-
length utterance-level descriptor. (2) The output descriptor should be compact
(i.e. low-dimensional), requiring little memory, to facilitate efficient
storage and retrieval. (3) The output descriptor should also be
discriminative, such that the distance between descriptors of different
speakers is larger than those of the same speaker.

To satisfy all the aforementioned properties, we use a modified ResNet in a
fully convolutional way to encode input 2D spectrograms, followed by a
NetVLAD/GhostVLAD layer for feature aggregation along the temporal axis. This
produces a fixed-length output descriptor. Intuitively, the VLAD layer can be
thought of as trainable discriminative clustering: every frame-level
descriptor will be softly assigned to different clusters, and residuals are
encoded as the output features. To allow efficient verification (i.e. low
memory, fast similarity computation), we further add a fully connected layer
for dimensionality reduction. Discriminative representations emerge because
the entire network is trained end-to-end for speaker identification. The
network is shown in Figure 1 and described in more detail in the following
paragraphs.

Fig. 1: Network architecture. It consists of two parts: _feature extraction_ ,
where a shared CNN is used to encode the spectrogram and extract frame-level
features, and _aggregation_ , which aggregates all the local descriptors into
a single compact representation of arbitrary length.

Module |  Input Spectrogram ($257\times T\times 1$) |  Output Size
---|---|---
_Thin ResNet_ |  conv2d, $7\times 7$, $64$ |  $257\times T\times 64$
max pool, $2\times 2$, stride ($2,2$) |  $128\times T/2\times 64$
$\begin{bmatrix}{\rm conv},1\times 1,48\\\ {\rm conv},3\times 3,48\\\ {\rm conv},1\times 1,96\end{bmatrix}\times 2$ |  $128\times T/2\times 96$
$\begin{bmatrix}{\rm conv},1\times 1,96\\\ {\rm conv},3\times 3,96\\\ {\rm conv},1\times 1,128\end{bmatrix}\times 3$ |  $64\times T/4\times 128$
$\begin{bmatrix}{\rm conv},1\times 1,128\\\ {\rm conv},3\times 3,128\\\ {\rm conv},1\times 1,256\end{bmatrix}\times 3$ |  $32\times T/8\times 256$
$\begin{bmatrix}{\rm conv},1\times 1,256\\\ {\rm conv},3\times 3,256\\\ {\rm conv},1\times 1,512\end{bmatrix}\times 3$ |  $16\times T/16\times 512$
max pool, $3\times 1$, stride ($2,2$) |  $7\times T/32\times 512$
conv2d, $7\times 1$, $512$ |  $1\times T/32\times 512$

Table 1: The thin-ResNet used for frame level feature extraction. ReLU and
batch-norm layers are not shown. Each row specifies the # of convolutional
filters, their sizes, and the # filters. This architecture has only 3 million
parameters compared to the standard ResNet-34 (22 million).

Feature Extraction. The first stage involves feature extraction from input
spectrograms. While any network can be used in our learning framework, we opt
for a modified ResNet with $34$ layers. Compared to the standard ResNet used
before by [3], we cut down the number of channels in each residual block,
making it a _thin_ ResNet-34 (Table 1).

NetVLAD. The second part of the network uses NetVLAD [16] to aggregate frame-
level descriptors into a single utterance-level vector. Here we provide a
brief overview of NetVLAD (for full details please refer to [16]).

The thin ResNet maps the input spectrogram ($R^{257\times T\times 1}$) to
frame-level descriptors with size $R^{1\times T/32\times 512}$. The NetVLAD
layer then takes dense descriptors as input and produces a single $K\times D$
matrix $V$, where $K$ refers to the number of chosen cluster, and $D$ refers
to the dimensionality of each cluster. Concretely, the matrix of descriptors
$V$ is computed using the following equation:

$\displaystyle
V(k,j)=\sum_{t=1}^{T/32}\frac{e^{w_{k}x_{t}+b_{k}}}{\sum_{k^{\prime}=1}^{K}e^{w_{k}^{\prime}x_{t}+b_{k^{\prime}}}}(x_{t}(j)-c_{k}(j))$
(1)

where $\\{w_{k}\\},\\{b_{k}\\}$ and $\\{c_{k}\\}$ are trainable parameters,
with $k\in[1,2,...,K]$. The first term corresponds to the soft-assignment
weight of the input vector $x_{i}$ for cluster $k$, while the second term
computes the residual between the vector and the cluster centre. The final
output is obtained by performing $L2$ normalisation and concatenation. To keep
computational and memory requirements low, dimensionality reduction is
performed via a Fully Connected (FC) layer, where we pick the output
dimensionality to be $512$. We also experiment with the recently proposed
GhostVLAD [17] layer, where some of the clusters are not included in the final
concatenation, and so do not contribute to the final representation, these are
referred to as ‘ghost clusters’ (we used _two_ in our implementation).
Therefore, while aggregating the frame-level features, the contribution of the
noisy and undesirable sections of a speech segment to normal VLAD clusters is
effectively down-weighted, as most of their weights have been assigned to the
‘ghost cluster’. For further details, please see [17].

##  3 Experiments

###  3.1 Datasets

We train our model end-to-end on the VoxCeleb2 [3] dataset (only on the ‘dev’
partition, this contains speech from 5,994 speakers) for identification and
test on the VoxCeleb1 verification test sets [3]. Note that the development
set of VoxCeleb2 is completely disjoint from the VoxCeleb1 dataset (i.e. no
speakers in common).

###  3.2 Training Loss

Besides the standard softmax loss, we also experiment with the additive margin
softmax (AM-Softmax) classification loss [25] during training. This loss is
known to improve verification performance by introducing a margin in the
angular space. The loss is given by the following equation:

$\displaystyle
L_{i}=-log\frac{e^{s(\cos\theta_{y_{i}}-m)}}{e^{s(\cos\theta_{y_{i}}-m)}+\sum_{j\neq{y_{i}}}e^{s\cos(\theta_{j})}}$
(2)

where $L_{i}$ refers to cost of assigning the sample to the correct class,
$\theta_{y}=arccos(w^{T}x)$ refers to the angle between sample features ($x$)
and the decision hyperplane ($w$), as both vectors have been L2 normalised.
The goal is therefore to minimise this angle by making $cos(\theta_{y_{i}})-m$
as large as possible, where $m$ refers to the angular margin. The hyper-
parameter $s$ controls the “temperature” of the softmax loss, producing higher
gradients to the well-separated samples (and further shrinking the intra-class
variance). We used the default values $m=0.4$ and $s=30$ [25].

###  3.3 Training Details

During training, we use a fixed size spectrogram corresponding to a 2.5 second
temporal segment, extracted randomly from each utterance. Spectrograms are
generated in a sliding window fashion using a hamming window of width 25ms and
step 10ms. We use a 512 point FFT, giving us 256 frequency components, which
together with the DC component of each frame gives a short-time Fourier
transform (STFT) of size $257\times 250$ (frequency $\times$temporal) out of
every 2.5 second crop. The spectrogram is normalised by subtracting the mean
and dividing by the standard deviation of all frequency components in a single
time step. No voice activity detection (VAD), or automatic silence removal is
applied. We use the Adam optimiser with an initial learning rate of $1e-3$,
and decrease the learning rate by 10 after every 36 epochs until convergence.

##  4 Results

In this section we first compare the performance of our NetVLAD and GhostVLAD
architectures trained using different losses to the state of the art, and then
investigate how performance varies with utterance length.

###  4.1 Verification on VoxCeleb1

| Front-end model | Loss | Dims | Aggregation | Training set | EER (%)
---|---|---|---|---|---|---
VoxCeleb1 test set |  |  |  |  |  |
Nagrani et al. [2] | I-vectors + PLDA | – | – | – | VoxCeleb1 | 8.8
Nagrani et al. [2] | VGG-M | Softmax | 1024 | TAP | VoxCeleb1 | 10.2
Cai et al. [15] | ResNet-34 | A-Softmax + PLDA | 128 | TAP | VoxCeleb1 | 4.46
Cai et al. [15] | ResNet-34 | A-Softmax + PLDA | 128 | SAP | VoxCeleb1 | 4.40
Cai et al. [15] | ResNet-34 | A-Softmax + PLDA | 128 | LDE | VoxCeleb1 | 4.48
Okabe et al. [13] | TDNN (x-vector) | Softmax | 1500 | TAP | VoxCeleb1 | 4.70
Okabe et al. [13] | TDNN (x-vector) | Softmax | 1500 | SAP | VoxCeleb1 | 4.19
Okabe et al. [13] | TDNN (x-vector) | Softmax | 1500 | ASP | VoxCeleb1 | 3.85
Hajibabaei et al. [19] | ResNet-20 | A-Softmax | 128 | TAP | VoxCeleb1 | 4.40
Hajibabaei et al. [19] | ResNet-20 | AM-Softmax | 128 | TAP | VoxCeleb1 | 4.30
Chung et al. [3] | ResNet-34 | Softmax + Contrastive | 512 | TAP | VoxCeleb2 | 5.04
Chung et al. [3] | ResNet-50 | Softmax + Contrastive | 512 | TAP | VoxCeleb2 | 4.19
Ours | Thin ResNet-34 | Softmax | 512 | TAP | VoxCeleb2 | 10.48
Ours | Thin ResNet-34 | Softmax | 512 | NetVLAD | VoxCeleb2 | 3.57
Ours | Thin ResNet-34 | AM-Softmax | 512 | NetVLAD | VoxCeleb2 | 3.32
Ours | Thin ResNet-34 | Softmax | 512 | GhostVLAD | VoxCeleb2 | 3.22
Ours | Thin ResNet-34 | AM-Softmax | 512 | GhostVLAD | VoxCeleb2 | 3.23
Ours (cleaned $\dagger$) | Thin ResNet-34 | Softmax | 512 | GhostVLAD | VoxCeleb2 | 3.24
VoxCeleb1-E |  |  |  |  |  |
Chung et al. [3] | ResNet-50 | Softmax + Contrastive | 512 | TAP | VoxCeleb2 | 4.42
Ours | Thin ResNet-34 | Softmax | 512 | GhostVLAD | VoxCeleb2 | 3.24
Ours (cleaned $\dagger$) | Thin ResNet-34 | Softmax | 512 | GhostVLAD | VoxCeleb2 | 3.13
VoxCeleb1-H |  |  |  |  |  |
Chung et al. [3] | ResNet-50 | Softmax + Contrastive | 512 | TAP | VoxCeleb2 | 7.33
Ours | Thin ResNet-34 | Softmax | 512 | GhostVLAD | VoxCeleb2 | 5.17
Ours (cleaned $\dagger$) | Thin ResNet-34 | Softmax | 512 | GhostVLAD | VoxCeleb2 | 5.06
Table 2: Results for verification on the original VoxCeleb1 test set [2] and
the extended and hard test sets (VoxCeleb-E and VoxCeleb-H) [3]. [15, 13, 19]
do not report results on the VoxCeleb-E and VoxCeleb-H) [3] test sets. TAP:
Temporal Average Pooling. SAP: Self-attentive Pooling Layer [15], $\dagger$
Cleaned up versions of the test lists have been released publically. We
encourage other researchers to evaluate on these lists.

The trained network is evaluated on three different test lists from the
VoxCeleb1 dataset: (1) the original VoxCeleb1 test list with 40 speakers; (2)
the extended VoxCeleb1-E list that uses the entire VoxCeleb1 (train and test
splits) for evaluation; and (3) the challenging VoxCeleb1-H list where the
test pairs are drawn from identities with the same gender and nationality. In
addition, we find that there are a small number of errors in the VoxCeleb1-E
and VoxCeleb1-H lists, and hence we evaluate on a cleaned up version of both
lists as well, which we release publically. The network is tested on the full
length of the test segment. We do not use any test time augmentation, which
could potentially lead to slight performance gains.

Table 2 compares the performance of our models to the current state-of-the-art
on the original VoxCeleb1 test set. Our model outperforms all previous
methods. With standard softmax loss and a NetVLAD aggregation layer, it
outperforms the original ResNet-based architecture [3] by a significant margin
(EER of 3.57% vs 4.19%) whilst requiring far fewer parameters (10 vs 26
million). By replacing the standard softmax with the additive margin softmax
(AM-Softmax), a further performance gain is achieved (3.32% EER). The
GhostVLAD layer, which excludes irrelevant information from the aggregation,
additionally makes a modest contribution to performance (3.22% EER). On the
challenging VoxCeleb1-H test set, we outperform the original ResNet-based
architecture [3] (EER of 5.17% vs 7.33%), which is by a larger margin than on
the original VoxCeleb1 test set. The most similar architecture to ours is the
dictionary based method of Cai et al.  [15], which we also outperform (EER of
3.22% vs 4.48% ). We note that training a softmax loss based on features from
temporal average pooling (TAP) yields extremely poor results (EER of 10.48%).
We conjecture that the features from TAP are typically good at optimizing the
inter-class difference (i.e., separating different speakers), but not good at
reducing the intra-class variation (i.e. making features of the same speaker
compact). Therefore, contrastive loss with online hard sample mining leads to
a significant performance boost, as demonstrated in [3] for TAP. It is
possible that this would also give a performance boost for NetVLAD/GhostVLAD
pooling.

###  4.2 Additional experiment on GhostVLAD

In Table 3, we study the effect of the number of clusters in the GhostVLAD
layer. Despite small differences in performance, we show that VLAD aggregation
is robust to the number of clusters and to two different loss functions.

Loss | Aggregation | Clusters | G_Clusters | EER (%)
---|---|---|---|---
Softmax | GhostVLAD | 8 | 2 | 3.22
AMSoftmax | GhostVLAD | 8 | 2 | 3.23
Softmax | GhostVLAD | 10 | 2 | 3.37
AMSoftmax | GhostVLAD | 10 | 2 | 3.34
Softmax | GhostVLAD | 12 | 2 | 3.30
Softmax | GhostVLAD | 14 | 2 | 3.31

Table 3: Results for verification on the original VoxCeleb1 test set. All
models used the same architecture (Thin ResNet-34), and we vary the number of
VLAD clusters and the loss function.

###  4.3 Probing verification based on length

Table 4 shows the effect of the length of the test segment on speaker
recognition performance. In order to provide a fair comparison on lengths up
to 6 seconds, we restricted the testing dataset (VoxCeleb1) to speech segments
that were 6 seconds or longer (87,010 segments or 56.7% of the total dataset).
To generate verification pairs, for each speaker in the VoxCeleb1 dataset
(1251 speakers in total), we randomly sample 100 positive pairs and 100
negative pairs, resulting in 25,020 verification pairs. During testing,
segments of length 2s, 3s, 4s, 5s and 6s are randomly cropped from each
verification pair. We repeat this process three times, and compute mean and
standard deviation.

Segment |  |  |  |  |
---|---|---|---|---|---
length(s) | 2 | 3 | 4 | 5 | 6
EER |  7.97$\pm$0.06 |  5.73$\pm$0.04 |  4.70$\pm$0.02 |  4.10$\pm$0.02 |  3.39$\pm$0.02
Table 4: Utterance length (in seconds) on performance.

As shown in Table 4, there is indeed a strong correlation between verification
performance and sequence length. For ‘in the wild’ sequences, some of the data
may be noise, silence, or speech from other speakers, and a single short
sequence may be unlucky and consist of predominantly these irrelevant signals.
As the temporal length increases, there is a higher chance of capturing
relevant speech signals from the actual speaker.

##  5 Conclusion

In this paper, we have proposed a powerful speaker recognition network, using
a ‘thin-ResNet’ trunk architecture, and a dictionary-based NetVLAD and
GhostVLAD layers to aggregate features across time that can be trained end-to-
end. The network achieves state-of-the-art performance on the popular
VoxCeleb1 test set for speaker recognition, whilst requiring fewer parameters
than previous methods. We have also shown the effect of utterance length on
performance, and concluded that for ‘in the wild’ data, a longer length is
beneficial.

Acknowledgements. Funding for this research is provided by the EPSRC Programme
Grant Seebibyte EP/M013774/1. AN is supported by a Google PhD Fellowship in
Machine Perception, Speech Technology and Computer Vision.

## References

  * [1] M. McLaren, L. Ferrer, D. Castan, and A. Lawson,  “The speakers in the wild (SITW) speaker recognition database,”  in INTERSPEECH, 2016.
  * [2] A. Nagrani, J. S. Chung, and A. Zisserman,  “VoxCeleb: a large-scale speaker identification dataset,”  in INTERSPEECH, 2017.
  * [3] J. S. Chung, A. Nagrani, and A. Zisserman,  “VoxCeleb2: Deep speaker recognition,”  in INTERSPEECH, 2018.
  * [4] M. Abadi, A. Agarwal, P. Barham, E. Brevdo, Z. Chen, C. Citro, G.S. Corrado, A. Davis, J. Dean, M. Devin, et al.,  “Tensorflow: Large-scale machine learning on heterogeneous distributed systems,”  arXiv preprint arXiv:1603.04467, 2016.
  * [5] A. Paszke, S. Gross, S. Chintala, G. Chanan, E. Yang, Z. DeVito, Z. Lin, A. Desmaison, L. Antiga, and A. Lerer,  “Automatic differentiation in pytorch,”  2017\.
  * [6] A. Vedaldi and K. Lenc,  “Matconvnet: Convolutional neural networks for matlab,”  in Proc. ACMM, 2015.
  * [7] C. Li, X. Ma, B. Jiang, X. Li, X. Zhang, X. Liu, Y. Cao, A. Kannan, and Z. Zhu,  “Deep speaker: an end-to-end neural speaker embedding system,”  arXiv preprint arXiv:1705.02304, 2017.
  * [8] L. Wan, Q. Wang, A. Papir, and I.L. Moreno,  “Generalized end-to-end loss for speaker verification,”  in 2018 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2018, pp. 4879–4883.
  * [9] Y. Lukic, C. Vogt, O. Dürr, and T. Stadelmann,  “Speaker identification and clustering using convolutional neural networks,”  in IEEE 26th International Workshop on Machine Learning for Signal Processing (MLSP). IEEE, 2016, pp. 1–6.
  * [10] I. Lopez-Moreno, J. Gonzalez-Dominguez, O. Plchot, D. Martinez, J. Gonzalez-Rodriguez, and P. Moreno,  “Automatic language identification using deep neural networks,”  in Acoustics, Speech and Signal Processing (ICASSP), 2014 IEEE International Conference on. IEEE, 2014, pp. 5337–5341.
  * [11] D. Snyder, D. Garcia-Romero, D. Povey, and S. Khudanpur,  “Deep neural network embeddings for text-independent speaker verification,”  Proc. Interspeech 2017, pp. 999–1003, 2017.
  * [12] S. Shon, H. Tang, and J. Glass,  “Frame-level speaker embeddings for text-independent speaker recognition and analysis of end-to-end model,”  arXiv preprint arXiv:1809.04437, 2018.
  * [13] K. Okabe, T. Koshinaka, and K. Shinoda,  “Attentive statistics pooling for deep speaker embedding,”  arXiv preprint arXiv:1803.10963, 2018.
  * [14] W. Cai, Z. Cai, X. Zhang, X. Wang, and M. Li,  “A novel learnable dictionary encoding layer for end-to-end language identification,”  arXiv preprint arXiv:1804.00385, 2018.
  * [15] W. Cai, J. Chen, and M. Li,  “Exploring the encoding layer and loss function in end-to-end speaker and language recognition system,”  arXiv preprint arXiv:1804.05160, 2018.
  * [16] R. Arandjelović, P. Gronat, A. Torii, T. Pajdla, and J. Sivic,  “NetVLAD: CNN architecture for weakly supervised place recognition,”  in Proc. CVPR, 2016.
  * [17] Y. Zhong, R. Arandjelović, and A. Zisserman,  “GhostVLAD for set-based face recognition,”  in Asian Conference on Computer Vision, ACCV, 2018.
  * [18] W. Cai, J. Chen, and M. Li,  “Analysis of length normalization in end-to-end speaker verification system,”  arXiv preprint arXiv:1806.03209, 2018.
  * [19] M. Hajibabaei and D. Dai,  “Unified hypersphere embedding for speaker recognition,”  arXiv preprint arXiv:1807.08312, 2018.
  * [20] G. Bhattacharya, J. Alam, and P. Kenny,  “Deep speaker embeddings for short-duration speaker verification,”  in Proc. Interspeech, 2017, pp. 1517–1521.
  * [21] D. Snyder, D. Garcia-Romero, G. Sell, D. Povey, and S. Khudanpur,  “X-vectors: Robust dnn embeddings for speaker recognition,”  ICASSP, Calgary, 2018.
  * [22] FA Chowdhury, Quan Wang, Ignacio Lopez Moreno, and Li Wan,  “Attention-based models for text-dependent speaker verification,”  arXiv preprint arXiv:1710.10470, 2017.
  * [23] Jinkun Chen, Weicheng Cai, Danwei Cai, Zexin Cai, Haibin Zhong, and Ming Li,  “End-to-end language identification using netfv and netvlad,”  in Proc. ISCSLP, 2018.
  * [24] D. Chen, S. Tsai, V. Chandrasekhar, G. Takacs, H. Chen, R. Vedantham, R. Grzeszczuk, and B. Girod,  “Residual enhanced visual vectors for on-device image matching,”  in Asilomar, 2011.
  * [25] F. Wang, W. Liu, H. Liu, and J. Cheng,  “Additive margin softmax for face verification,”  arXiv preprint arXiv:1801.05599, 2018.


* Xie  _et al._ (1999) M. H. Xie, S. M. Seutter, W. K. Zhu, L. X. Zheng, Huasheng Wu, and S. Y. Tong, “Anisotropic step-flow growth and island growth of GaN(0001) by molecular beam epitaxy,” Phys. Rev. Lett. 82, 2749–2752 (1999).
  * Giesen (2001) Margret Giesen, “Step and island dynamics at solid/vacuum and solid/liquid interfaces,” Prog. Surf. Sci. 68, 1–154 (2001).
  * Verma (1951) Ajit Ram Verma, “CI. Observations on carborundum of growth spirals originating from screw dislocations,” Philos. Mag. 42, 1005–1013 (1951).
  * Sunagawa and Bennema (1979) I. Sunagawa and P. Bennema, “Modes of vibrations in step trains: Rhythmical bunching,” J. Cryst. Growth 46, 451–457 (1979).
  * van der Hoek  _et al._ (1982) B. van der Hoek, J. P. van der Eerden, and K. Tsukamoto, “Interpretation of double spirals on silicon carbide,” J. Cryst. Growth 58, 545–553 (1982).
  * Heying  _et al._ (1999) B. Heying, E. J. Tarsa, C. R. Elsass, P. Fini, S. P. DenBaars, and J. S. Speck, “Dislocation mediated surface morphology of GaN,” J. Appl. Phys. 85, 6470–6476 (1999).
  * S Vézian  _et al._ (2001) S. S Vézian, J. Massies, F. Semond, and N. Grandjean, “Surface morphology of GaN grown by molecular beam epitaxy,” Mater. Sci. Eng. B 82, 56–58 (2001).
  * Zauner  _et al._ (2002) A. R. A Zauner, E. Aret, W. J. P. van Enckevort, J. L. Weyher, S. Porowski, and J. J. Schermer, “Homo-epitaxial growth on the N-face of GaN single crystals: the influence of the misorientation on the surface morphology,” J. Cryst. Growth 240, 14–21 (2002).
  * Xie  _et al._ (2006) M. H. Xie, M. Gong, E. K. Y. Pang, H. S. Wu, and S. Y. Tong, “Origin of triangular island shape and double-step bunching during GaN growth by molecular-beam epitaxy under excess Ga conditions,” Phys. Rev. B 74, 085314 (2006).
  * Krukowski  _et al._ (2007) S. Krukowski, P. Kempisty, P. Str ak, G. Nowak, R. Czernecki, M. Leszczynski, T. Suski, M. Bockowski, and I. Grzegory, “Modelling the growth of nitrides in ammonia-rich environment,” Cryst. Res. Technol. 42, 1281–1290 (2007).
  * Zheng  _et al._ (2008) Hao Zheng, M. H. Xie, H. S. Wu, and Q. K. Xue, “Kinetic energy barriers on the GaN(0001) surface: A nucleation study by scanning tunneling microscopy,” Phys. Rev. B 77, 045303 (2008).
  * Turski  _et al._ (2013) H. Turski, M. Siekacz, Z. R. Wasilewski, M. Sawicka, S. Porowski, and C. Skierbiszewski, “Nonequivalent atomic step edges - role of gallium and nitrogen atoms in the growth of InGaN layers,” J. Cryst. Growth 367, 115–121 (2013).
  * Lin  _et al._ (2013) Chia-Hung Lin, Tetsuya Akasaka, and Hideki Yamamoto, “Nucleus and spiral growth of N-face GaN(000$\overline{1}$) obtained by selective-area metalorganic vapor phase epitaxy,” Appl. Phys Express 6, 035503 (2013).
  * Pristovsek  _et al._ (2017) Markus Pristovsek, Konrad Bellman, Frank Mehnke, Joachim Stellmach, Tim Wernicke, and Michael Kneissl, “Surface reconstructions of (0001) AlN during metal-organic vapor phase epitaxy,” Phys. Status Solidi B 254, 1600711 (2017).
  * Chen  _et al._ (2002) Yefan Chen, Hang-Ju Ko, Soon-Ku Hong, Takafumi Yao, and Yusaburo Segawa, “Morphology evolution of ZnO(000$\overline{1}$) surface during plasma-assisted molecular-beam epitaxy,” Appl. Phys. Lett. 80, 1358–1360 (2002).
  * Załuska-Kotur  _et al._ (2011) Magdalena A. Załuska-Kotur, Filip Krzyżewski, and Stanisław Krukowski, “Double step structure and meandering due to the many body interaction at GaN(0001) surface in N-rich conditions,” J. Appl. Phys. 109, 023515 (2011).
  * Załuska-Kotur  _et al._ (2010) Magdalena A. Załuska-Kotur, Filip Krzyżewski, and Stanisław Krukowski, “Surface patterns due to step flow anisotropy formed in crystal growth process,” J. Non-Cryst. Solids 356, 1935–1939 (2010).
  * Xu  _et al._ (2017) Dongwei Xu, Peter Zapol, G. Brian Stephenson, and Carol Thompson, “Kinetic monte carlo simulations of GaN homoepitaxy on c- and m-plane surfaces,” J. Chem. Phys. 146, 144702 (2017).
  * Chugh and Ranganathan (2017a) Manjusha Chugh and Madhav Ranganathan, “Lattice kinetic monte carlo simulation study of the early stages of epitaxial GaN(0001) growth,” Appl. Surf. Sci. 422, 1120–1128 (2017a).
  * Akiyama  _et al._ (2020a) Toru Akiyama, Takumi Ohka, Kohji Nakamura, and Tomonori Ito, “Ab initio study for adsorption and desorption behavior at step edges of GaN(0001) surface,” J. Cryst. Growth 532, 125410 (2020a).
  * Akiyama  _et al._ (2020b) Toru Akiyama, Takumi Ohka, Kohji Nakamura, and Tomonori Ito, “Ab initio study for adsorption and desorption behavior at step edges of AlN(0001) and GaN(0001) surfaces,” Jpn. J. Appl. Phys. 59, SGGK03 (2020b).
  * Avery  _et al._ (1997) A. R. Avery, H. T. Dobbs, D. M. Holmes, B. A. Joyce, and D. D. Vvedensky, “Nucleation and growth of islands on GaAs surfaces,” Phys. Rev. Lett. 79, 3938–3941 (1997).
  * Jo  _et al._ (2012) Masafumi Jo, Takaaki Mano, Marco Abbarchi, Takashi Kuroda, Yoshiki Sakuma, and Kasuaki Sakoda, “Self-limiting growth of hexagonal and triangular quantum dots on (111)A,” Cryst. Growth Des. 12, 1411–1415 (2012).
  * Kalff  _et al._ (1998) Matthias Kalff, George Comsa, and Thomas Michely, “How sensitive is epitaxial growth to adsorbates?” Phys. Rev. Lett. 81, 1255–1258 (1998).
  * Yin  _et al._ (2009) Cong Yin, Xi-Jing Ning, Jun Zhuang, Yi-Qun Xie, Xiu-Fang Gong, Xiang-Xi Ye, Chen Ming, and Yun-Fei Jin, “Shape prediction of two-dimensional adatom islands on crystal surfaces during homoepitaxial growth,” Appl. Phys. Lett. 94, 183107 (2009).
  * Robinson (1986) Ian K. Robinson, “Crystal truncation rods and surface roughness,” Phys. Rev. B 33, 3830 (1986).
  * Walkosz  _et al._ (2012) Weronika Walkosz, Peter Zapol, and G. Brian Stephenson, “Metallicity of InN and GaN surfaces exposed to NH3,” Phys. Rev. B 85, 033308 (2012).
  * Otte and Crocker (1965) H. M. Otte and A. G. Crocker, “Crystallographic formulae for hexagonal lattices,” Phys. Status Solidi 9, 441–450 (1965).
  * Munkholm and Brennan (1999) A. Munkholm and S. Brennan, “Influence of miscut on crystal truncation rod scattering,” J. Appl. Crystallogr. 32, 143–153 (1999).
  * Trainor  _et al._ (2002) Thomas P. Trainor, Peter J. Eng, and Ian K. Robinson, “Calculation of crystal truncation rod structure factors for arbitrary rational surface terminations,” J. Appl. Crystallogr. 35, 696–701 (2002).
  * Waasmaier and Kirfel (1995) D. Waasmaier and A. Kirfel, “New analytical scattering-factor functions for free atoms and ions,” Acta Crystallogr. Sec. A 51, 416–431 (1995).
  * Henke  _et al._ (1993) B. L. Henke, E. M. Gullikson, and J. C. Davis, “X-ray interactions: photoabsorption, scattering, transmission, and reflection at E=50-30000 eV, Z=1-92,” At. Data Nucl. Data Tables 54, 181 (1993), updates at https://henke.lbl.gov/optical_constants/asf.html.
  * Ju  _et al._ (2017) Guangxu Ju, Matthew J. Highland, Angel Yanguas-Gil, Carol Thompson, Jeffrey A. Eastman, Hua Zhou, Sean M. Brennan, G. Brian Stephenson, and Paul H. Fuoss, “An instrument for in situ coherent x-ray studies of metal-organic vapor phase epitaxy of III-nitrides,” Rev. Sci. Instrum. 88, 035113 (2017).
  * Ju  _et al._ (2018) Guangxu Ju, Matthew J. Highland, Carol Thompson, Jeffrey A. Eastman, Paul H. Fuoss, Hua Zhou, Roger Dejus, and G. Brian Stephenson, “Characterization of the x-ray coherence properties of an undulator beamline at the Advanced Photon Source,” J. of Synchrotron Radiat. 25, 1036–1047 (2018).
  * Ju  _et al._ (2019) Guangxu Ju, Dongwei Xu, Matthew J. Highland, Carol Thompson, Hua Zhou, Jeffrey A. Eastman, Paul H. Fuoss, Peter Zapol, Hyunjung Kim, and G. Brian Stephenson, “Coherent x-ray spectroscopy reveals the persistence of island arrangements during layer-by-layer growth,” Nat. Phys. 15, 589–594 (2019).
  * Note (1) GANKIBANTM from SixPoint Materials, Inc., spmaterials.com.
  * Murty  _et al._ (2000) M. V. Ramana Murty, P. Fini, G. B. Stephenson, Carol Thompson, J. A. Eastman, A. Munkholm, O. Auciello, R. Jothilingam, S. P. DenBaars, and J. S. Speck, “Step bunching on the vicinal GaN (0001) surface,” Phys. Rev. B 62, R10661 (2000).
  * Kempisty and Kangawa (2019) Pawel Kempisty and Yoshihiro Kangawa, “Evolution of the free energy of the GaN (0001) surface based on first-principles phonon calculations,” Phys. Rev. B 100, 085304 (2019).
  * Guin  _et al._ (2020) L. Guin, M. E. Jabbour, L. Shaabani-Ardali, L. Benoit-Maréchal, and N. Triantafyllidis, “Stability of vicinal surfaces: Beyond the quasistatic approximation,” Phys. Rev. Lett. 124, 036101 (2020).
  * Li  _et al._ (2016) Yuan Li, Xuejiang Chen, and Juan Su, “Study on formation of step bunching on 6H-SiC (0001) surface by kinetic Monte Carlo method,” Appl. Surf. Sci. 371, 242–247 (2016).
  * Bellmann  _et al._ (2017) Konrad Bellmann, Udo W. Pohl, Christian Kuhn, Tim Wernicke, and Michael Kneissl, “Controlling the morphology transition between step-flow growth and step-bunching growth,” J. Cryst. Growth 478, 187–192 (2017).
  * Dufay  _et al._ (2007) Matthieu Dufay, Thomas Frisch, and Jean-Marc Debierre, “Role of step-flow advection during electromigration-induced step bunching,” Phys. Rev. B 75, 241304 (2007).
  * Pierre-Louis (2003a) O. Pierre-Louis, “Step bunching with general step kinetics: stability analysis and macroscopic models,” Surf. Sci. 529, 114–134 (2003a).
  * Pimpinelli and Videcoq (2000) A. Pimpinelli and A. Videcoq, “Novel mechanism for the onset of morphological instabilities during chemical vapour epitaxial growth,” Surf. Sci. 445, L23–L28 (2000).
  * Pierre-Louis and Métois (2004) O. Pierre-Louis and J.-J. Métois, “Kinetic step pairing,” Phys. Rev. Lett. 93, 165901 (2004).
  * Patrone  _et al._ (2010) Paul N. Patrone, T. L. Einstein, and Dionisios Margetis, “One-dimensional model of interacting-step fluctuations on vicinal surfaces: Analytical formulas and kinetic Monte Carlo simulations,” Phys. Rev. E 82, 061601 (2010).
  * Ranguelov  _et al._ (2007) Bogdan Ranguelov, M. S. Altman, and Ivan Markov, “Critical terrace width for step flow growth: Effect of attachment-detachment asymmetry and step permeability,” Phys. Rev. B 75, 245419 (2007).
  * Hanada (2019) Takashi Hanada, “Thermodynamic model for metalorganic vapor-phase epitaxy of N-polar group-III nitrides in step-flow growth mode: Hydrogen, competitive adsorption, and configuration entropy,” Phys. Rev. Materials 3, 103404 (2019).
  * Pierre-Louis (2003b) O. Pierre-Louis, “Phase field models for step flow,” Phys. Rev. E 68, 021604 (2003b).
  * Reeber and Wang (2000) Robert R. Reeber and Kai Wang, “Lattice parameters and thermal expansion of GaN,” J. Mater. Res. 15, 40–44 (2000).
  * Zywietz  _et al._ (1998) Tosja Zywietz, Jörg Neugebauer, and Matthias Scheffler, “Adatom diffusion at GaN(0001 and (000$\overline{1}$) surfaces,” Appl. Phys. Lett. 73, 487–489 (1998).
  * Chugh and Ranganathan (2017b) Manjusha Chugh and Madhav Ranganathan, “Adsorbate interactions on the GaN(0001) surface and their effect on diffusion barriers and growth morphology,” Phys. Chem. Chem. Phys. 19, 2111–2123 (2017b).
  * González-Hernández  _et al._ (2011) Rafael González-Hernández, William López-Pérez, María Guadalupe Moreno-Armenta, and Jairo Arbey Rodríguez M, “Adsorption and diffusion of 3d transition metal atoms on the GaN(0001) surface,” J. Appl. Phys. 110, 083712 (2011).
  * Koleske  _et al._ (2014) D. D. Koleske, S. R. Lee, M. H. Crawford, K. C. Cross, M. E. Coltrin, and J. M. Kempisty, “Connection between GaN and InGaN growth mechanisms and surface morphology,” J. Cryst. Growth 391, 85–96 (2014).
  * Shewmon (1989) Paul Shewmon,  _Diffusion in Solids_, 2nd ed. (Springer, 1989).
  * Van de Walle and Neugebauer (2002a) Chris G. Van de Walle and J. Neugebauer, “First-principles surface phase diagram for hydrogen on GaN surfaces,” Phys. Rev. Lett. 88, 066103 (2002a).
  * Van de Walle and Neugebauer (2002b) Chris G. Van de Walle and J. Neugebauer, “Role of hydrogen in surface reconstructions and growth of GaN,” J. Vac. Sci. Technol. B 20, 1640–1646 (2002b).
  * Chase (1998) Malcolm W. Chase, Jr., “NIST-JANAF thermochemical tables (4th edition),” J. Phys. Chem. Ref. Data, Monograph 9 , 1–1151 (1998).
  * Ambacher  _et al._ (1996) O. Ambacher, M. S. Brandt, R. Dimitrov, T. Metzger, M. Stutzmann, R. A. Fischer, A. Miehr, A. Bergmaier, and G. Dollinger, “Thermal stability and desorption of Group III nitrides prepared by metal organic chemical vapor deposition,” J. Vac. Sci. Technol. B 14, 3532–3542 (1996).
  * Jiang  _et al._ (2008) Fan Jiang, A. Munkholm, R.-V. Wang, S. K. Streiffer, Carol Thompson, P. H. Fuoss, K. Latifi, K. R. Elder, and G. B. Stephenson, “Spontaneous oscillations and waves during chemical vapor deposition of InN,” Phys. Rev. Lett. 101, 086102 (2008), (See online supplemental for trimethylindium (TMI) partial pressures required to condense InN and elemental In).
  * Stephenson  _et al._ (1999) G. B. Stephenson, J. A. Eastman, O. Auciello, A. Munkholm, Carol Thompson, P. H. Fuoss, P. Fini, S. P. DenBaars, and J. S. Speck, “Real-time x-ray scattering studies of surface structure during metalorganic chemical vapor deposition of GaN,” MRS Bull. 24[1], 21 (1999).
  * Stampfl and Van de Walle (1999) C. Stampfl and C. G. Van de Walle, “Density-functional calculations for III-V nitrides using the local-density approximation and the generalized gradient approximation,” Phys. Rev. B 59, 5521–5535 (1999).
  * Minikayev  _et al._ (2015) R. Minikayev, W. Paszkowicz, P. Piszora, M. Knapp, C. Bähtz, and S. Podsiadło, “Thermal expansion of polycrystalline gallium nitride: an x-ray diffraction study,” X-Ray Spectrom. 44, 382–388 (2015).
  * Perret  _et al._ (2014) Edith Perret, M. J. Highland, G. B. Stephenson, S. K. Streiffer, P. Zapol, P. H. Fuoss, A. Munkholm, and Carol Thompson, “Real-time x-ray studies of crystal growth modes during metal-organic vapor phase epitaxy of GaN on c- and m-plane single crystals,” Appl. Phys. Lett. 105, 051602 (2014).
  * Koleske  _et al._ (2005) D. D. Koleske, M. E. Coltrin, and M. J. Russell, “Using optical reflectance to measure GaN nucleation layer decomposition layer kinetics,” J. Cryst. Growth 279, 37–54 (2005).
  * Tapping and Reilly (1986) J. Tapping and M. L. Reilly, “Index of refraction of sapphire between 24 and 1060 C for wavelengths of 633 and 799 nm,” J. Opt. Soc. Am. A 3, 610–616 (1986).
  * Touloulian  _et al._ (1977) Y. S. Touloulian, R. K. Kirby, R. E. Taylor, and T. Y. Lee,  _Thermal Expansion: Nonmetallic Solids_, Vol. 13 (Springer, 1977) pp. 154–389.
  * Koleske  _et al._ (2001) D. D. Koleske, A. E. Wickenden, R. L. Henry, Culbertson. J. C., and M. E. Twigg, “GaN decomposition in H2 and N2 at MOVPE temperatures and pressures,” J. Cryst. Growth 223, 466–483 (2001).