# Speech Command Recognition in Computationally Constrained Environments with
a Quadratic Self-organized Operational Layer ††thanks: This work received
funding from the European Union’s Horizon 2020 research and innovation
programme under grant agreement No 871449 (OpenDR). This publication reflects
the authors’ views only. The European Commission is not responsible for any
use that may be made of the information it contains.

Mohammad Soltanian  Department of Computing Sciences
Tampere University
Tampere, Finland
<EMAIL_ADDRESS>Junaid Malik  Department of Computing Sciences
Tampere University
Tampere, Finland
<EMAIL_ADDRESS>Jenni Raitoharju  Programme for Environmental
Information
Finnish Environment Institute
Jyväskylä, Finland
<EMAIL_ADDRESS>Alexandros Iosifidis  Department of Engineering
Aarhus University
Aarhus, Denmark
<EMAIL_ADDRESS>Serkan Kiranyaz  Electrical Engineering Department
Qatar University
Doha, Qatar
<EMAIL_ADDRESS>Moncef Gabbouj  Department of Computing Sciences
Tampere University
Tampere, Finland
<EMAIL_ADDRESS>
###### Abstract

Automatic classification of speech commands has revolutionized human computer
interactions in robotic applications. However, employed recognition models
usually follow the methodology of deep learning with complicated networks
which are memory and energy hungry. So, there is a need to either squeeze
these complicated models or use more efficient light-weight models in order to
be able to implement the resulting classifiers on embedded devices. In this
paper, we pick the second approach and propose a network layer to enhance the
speech command recognition capability of a lightweight network and demonstrate
the result via experiments. The employed method borrows the ideas of Taylor
expansion and quadratic forms to construct a better representation of features
in both input and hidden layers. This richer representation results in
recognition accuracy improvement as shown by extensive experiments on Google
speech commands (GSC) and synthetic speech commands (SSC) datasets.

##  I Introduction

Automatic Speech recognition (ASR) is the art and science of automatic
identification of speech. Speech commands recognition (SCR) is a subcategory
of ASR in which the machine identifies short spoken commands. It is a key
technology which enables machines to understand human commands. Possible
applications include smart-phones with mobile assistants, robots which follow
human spoken instructions, and smart home assistants, like Amazon Echo and
Google Home. These kinds of applications are increasingly affecting human
machine interactions. There are also many use cases for this technology in
both medicine and education, e.g., for blind or handicapped people [1, 2]. As
many of these applications rely on small sized devices which have limited
computational capacities, fast and lightweight implementations are needed to
work smoothly in real time [3].

Classical SCR approaches rely on extracting discriminative features from the
raw speech. Frequency domain features like fast Fourier transform are among
the most widely used features [4, 5] to convert the time domain speech into
frequency space. After feature extraction, an acoustic model such as hidden
Markov model (HMM) was classically used to represent the sequence of words of
phonemes [6]. Combination of neural networks and HMM-based approaches have
also been in use since more than 20 years ago [7, 8].

Recently, deep learning-based approaches have dominated, in terms of
performance accuracy, over conventional machine learning methods for many
different applications [9, 10, 11]. Neural networks like long short term
memory (LSTM) have raised the accuracy of natural language processing systems
beyond the traditional approaches [12, 13]. Also, the neural network based
approaches have dominated the whole SCR history and have been accepted as the
state of the art tool by the speech recognition society [14, 15, 16, 17]. In
the meanwhile, several different neural network architectures like
convolutional neural networks (CNNs), LSTMs, and gated recurrent units (GRUs)
have been successfully applied to SCR [15, 16, 17]. Direct application of CNNs
on the raw audio have been explored in [18, 19, 20]. This is while the low
computational and memory capacities of the embedded devices are still serious
barriers to fully exploit the benefits of highly accurate complicated deep
networks.

To tackle the problem of computational complexity on embedded devices, there
are usually two general approaches. One can squeeze the state of the art deep
models to reduce the computational costs without severely degrading the
performance accuracy [21, 22]. Alternatively, more efficient deep networks and
methods can be employed from the beginning in order to achieve better accuracy
rates with lower computational costs [23, 24]. The second approach especially
becomes important when training examples are scarce or the input signals are
of low quality or resolution. Using either of these ways, one can implement a
reasonable classifier on energy constrained embedded devices.

There are many studies focusing on squeezing the complicated networks or
building lightweight networks from the beginning. Progressive operational
perceptrons (POPs)[25] try to generalize the conventional multilayer
perceptrons (MLPs) using generalized operational perceptrons (GOPs). GOPs have
a distinct set of operators, to mimic the synaptic connections of biological
neurons, which enables them to better learn the hidden semantics of data.
PyGOP [26] is a library for implementation of GOPs in both CPU and GPU.
Authors in [27] accelerate the learning process in POPs employing the concept
of memory and information from previously learned layers. The idea of
biologically inspired neurons is further explored in [28], where each neuron
in the network can have its distinct set of operators, which results in a
compact network.

In [29], a generic neural network layer structure employing a multi-linear
projection is proposed, which requires several times less memory compared to
traditional CNNs, while achieving better performance. Authors in [15] used
deep networks with dilated convolutions to achieve a high accuracy with lower
number of parameters. SincNet Kernels [19] at the first convolutional layer of
a CNN is another method to shrink the CNN by lowering the number of learnable
parameters. This method is used for speaker recognition in [18] and later for
phoneme recognition [19]. Depth-wise separable convolutions (DSCconvs) have
also been introduced firstly in image processing field [30, 31] and later in
SCR, keyword spotting [15], and machine translation [32] in order to reduce
the amount of required computations.

In this paper, we propose a new network layer to enhance the speech command
recognition capability using a lightweight network. To this end, we use a
novel combination of self-organized operational neural networks (SelfONNs)
[24] and quadratic forms kernels [33]. The experiments on Google speech
commands (GSC) [34] and synthetic speech commands (SSC)
datasets111https://kaggle.com/jbuchner/synthetic-speech-commands-dataset show
the effectiveness of the proposed approach for recognition of speech commands
especially in lightweight networks. To the best of our knowledge, we are the
first to propose a combination of SelfONNs and quadratic kernels to make a new
network layer type to boost the performance of a lightweight CNN. We are also
the first to apply concepts of SelfONN and quadratic kernels on the task of
SCR.

In the following, we describe the most relevant work in Section II. The
proposed quadratic SelfONN approach to the problem of SCR is introduced in
Section III, and the experiments are discussed in Section IV. Finally, the
paper is concluded in Section V.

##  II Related Work

###  II-A Self-organized Operational Neural Networks

Operational Neural Networks (ONNs) [23] are heterogeneous networks with a
generalized neuron model that includes a dictionary of nodal operators
(instead of just an ordinary convolution) to boost the performance of the
conventional CNNs. However, they need a greedy search to find the optimal
operators which result in the best accuracy. The performance is also highly
dependent on the selection of the initial operator set dictionary. In order to
resolve these issues, SelfONNs with generative neurons were recently proposed
in [24], which are able to adapt the nodal operators of each neuron during the
training phase. This removes the need for a fixed operator set and greedy
search within that library to find the optimal set of operators. The authors
show, via extensive experiments, the superiority of SelfONNs over ONNs in four
applications: image synthesis, image denoising, face segmentation, and image
transformation.

Let us denote by $\bm{X}\in\mathbb{R}^{H\times W\times C_{in}}$ the input
tensor to a layer and by $\bm{X}_{(i,j)}\in\mathbb{R}^{h\times w\times
C_{in}}$ a sub-tensor of it centered at the position $(i,j)$. Let us also
denote by $\bm{W}_{c_{out}}\in\mathbb{R}^{h\times w\times
C_{in}},\>c_{out}=1,\dots,C_{out}$ the $c_{out}$-th filter of the layer.
Convolutional neurons in CNNs convolve $\bm{X}$ with $\bm{W}_{c_{out}}$ and
add an offset, which corresponds to applying the following calculation for
every position $(i,j)$ of the input tensor $\bm{X}$:

$\displaystyle\bm{Y}_{c_{out}}(i,j)$
$\displaystyle=\sum_{k,m,c_{in}=1}^{h,w,C_{in}}\bm{W}_{c_{out}}(k,m,c_{in})\bm{X}_{(i,j)}(k,m,c_{in})+b_{c_{out}}$
$\displaystyle=\bm{w}_{c_{out}}^{T}\bm{x}_{(i,j)}+b_{c_{out}},$ (1)

where $\bm{Y}_{c_{out}}(i,j)$ is the $(i,j)$-th element of the output feature
map $\bm{Y}_{c_{out}}$, $b_{c_{out}}$ is the bias term, and $\bm{x}_{(i,j)}$
and $\bm{w}_{c_{out}}$ are the vectorized versions of $\bm{X}_{(i,j)}$ and
$\bm{W}_{c_{out}}$, respectively. The feature maps
$\bm{Y}_{c_{out}},\>c_{out}=1,\dots,C_{out}$ are concatenated to form the
tensor222We assume that appropriate zero-padding is used.
$\bm{Y}\in\mathbb{R}^{H\times W\times C_{out}}$, and an element-wise
activation function is applied to produce the output of the layer.

In ONNs, the transformation in Eq. (1) is generalized to the following:

$\bm{Y}_{c_{out}}(i,j)=\Psi(\bm{x}_{(i,j)},\bm{w}_{c_{out}})+b_{c_{out}},$ (2)

where $\Psi(\cdot)$ is an arbitrary (differentiable) nodal function and may be
a combination of different functions. The selection of $\Psi(\cdot)$ is
conducted by following a search strategy during training, which leads to a
slow training process. When $\Psi(\cdot)$ is determined to be the dot-product
of its arguments, the ONN layer is equal to an ordinary convolutional layer.

In SelfONNs, instead of searching for an optimal function $\Psi(\cdot)$ during
training, a suitable function is approximated with a truncated Taylor series
expansion:

$\displaystyle\Psi(\bm{x}_{(i,j)},\bm{w}_{c_{out},1},\dots,\bm{w}_{c_{out},Q})=\bm{w}_{c_{out},0}^{T}\bm{1}+\bm{w}_{c_{out},1}^{T}\bm{x}_{(i,j)}$
$\displaystyle+\bm{w}_{c_{out},2}^{T}\bm{x}_{(i,j)}^{2}+...+\bm{w}_{c_{out},Q}^{T}\bm{x}_{(i,j)}^{Q}=\sum_{q=0}^{Q}\bm{w}_{c_{out},q}^{T}\bm{x}_{(i,j)}^{q},$
(3)

where $\bm{1}$ is a vectors of ’1’s, in $\bm{x}_{(i,j)}^{q}$ the power is
applied in an element-wise manner and $\bm{w}_{c_{out},q}$ are the learnable
weights interacting with $\bm{x}_{(i,j)}^{q}$. Thus, each neuron of SelfONN
performs the following transformation:

$\bm{Y}_{c_{out}}(i,j)=\sum_{q=1}^{Q}\bm{w}_{c_{out},q}^{T}\bm{x}_{(i,j)}^{q}+b_{c_{out}},$
(4)

where the term $\bm{w}_{c_{out},0}$ is absorbed in the learnable bias
$b_{c_{out}}$. As the goal is not to estimate any particular function
$\Psi(\cdot)$, but to learn the most suitable nodal function for the learning
problem at hand, $\bm{w}_{c_{out},q},\>q=1,\dots,Q$ are learned via gradient-
based optimization.

###  II-B Quadratic-form Kernels

The core idea behind the network layer based on the quadratic form expansion
[33] is that it generalizes the linear convolution by taking into account the
cross correlation between the input elements within the receptive field of the
layer kernel, i.e.:

$\bm{Y}_{c_{out}}(i,j)=\bm{x}_{(i,j)}^{T}\bm{\Omega}_{c_{out}}\bm{x}_{(i,j)}+\bm{w}_{c_{out}}^{T}\bm{x}_{(i,j)}+b_{c_{out}}.$
(5)

Because the quadratic term is performed in a channel-wise manner,
$\bm{\Omega}_{c_{out}}$ is a block-diagonal matrix formed by upper-triangular
blocks. The diagonal block interacting with the elements of $\bm{X}$ in its
$c_{in}$-th channel, i.e. $\bm{\Omega}_{c_{out}}^{(c_{in},c_{in})}$ is a
Volterra kernel of the form:

$\bm{\Omega}_{c_{out}}^{(c_{in},c_{in})}=\begin{pmatrix}\omega_{1,1}&\omega_{1,2}&\cdots&\omega_{1,n}\\\
0&\omega_{2,2}&\cdots&\omega_{2,n}\\\ \vdots&\vdots&\ddots&\vdots\\\
0&0&\cdots&\omega_{n,n}\end{pmatrix},$ (6)

where $n=hw$.

##  III Proposed Quadratic Self-organized Operational Neural Network method

The concepts of SelfONN and quadratic convolution have not yet been applied to
SCR problem. We combine these two layer types into a new type of layer which
inherits advantages of both SelfONN (by a generalized convolution) and
quadratic kernels (by modeling cross correlation of input feature points). To
do so, we merge the Taylor series expansion in Eq. (3) with the quadratic-form
convolution in Eq. (5) as follows:

$\displaystyle\bm{Y}_{c_{out}}(i,j)$
$\displaystyle=\sum_{q=1}^{Q}(\bm{x}_{(i,j)}^{q})^{T}\bm{\Omega}_{c_{out},q}\bm{x}_{(i,j)}^{q}+\sum_{q=1}^{Q}\bm{w}_{c_{out},q}^{T}\bm{x}_{(i,j)}^{q}+b_{c_{out}}.$
(7)

where two sets of learnable weights $\bm{\Omega}_{c_{out},q}$ and
$\bm{w}_{c_{out},q}$ are interacting with $\bm{x}_{(i,j)}^{q}$. The matrix
$\bm{\Omega}_{c_{out},q}$ in Eq. (7) is again a block-diagonal matrix as in
Eq. (5) but, instead of using upper-triangular blocks, all values of its
diagonal blocks are learnable. The flow diagram of the proposed layer is shown
in Figure 1.

Figure 1: Flow diagram of the proposed quadratic-SelfONN based layer: The
spectrogram, MFCC, or any feature map is fed as the input. The input is then
expanded as a power series, and linear and quadratic convolutions are applied.
The results are summed together and the bias is added.

##  IV Experiments

###  IV-A Datasets

The proposed model is evaluated on Google Speech Commands (GSC) [34] and
Synthetic Speech Commands (SSC) 333https://kaggle.com/jbuchner/synthetic-
speech-commands-dataset datasets. GSC is the most widely used benchmark for
SCR systems. The training set of it consists of more than 60,000 audio files
categorized into 32 directories with the folder names being the label of the
audio clips. A subset of the dataset, which is commonly used, e.g., in Kaggle
competition, contains ”on”, ”off”, ”yes”, ”no”, ”left”, ”right”, ”up”, ”down”,
”stop”, and ”go” commands. This is the subset we use for our evaluations in
the current work. Each audio file in GSC is a voice clip with length of almost
1 second. Every audio file is in the form of 16-bit little-endian PCM-encoded
WAVE file at a sampling frequency of 16000 Hz.

SSC is automatically generated and is actually text-to-speech counterpart for
the GSC. It contains two versions, the normal version and very noisy version,
each of which consists of 41849 audio files categorized into 30 classes. Since
GSC is a low noise dataset, the very noisy version of SSC is employed in the
experiments to assess the performance of the proposed method on noisy data as
well. The same subset of classes as for GSC is used SSC in our experiments. In
contrary to GSC, the audio files in SSC are not segmented into different
train, test, and validation subsets. So, we randomly partition each category
into train (80 %), test (10 %), and validation (10 %).

The accuracy of the proposed method against the baseline is evaluated based on
the number of correctly predicted labels in the test set divided by the total
number of test samples, for both GSC and SSC datasets.

###  IV-B Data Pre-processing

The 1D raw input audio is converted to a 2D signal before being fed to the
network. Actually, each input speech wave is segmented to 30 milliseconds
windows with 10 milliseconds overlap. Then, Mel-frequency cepstral
coefficients (MFCCs) are extracted from each window. Since the configuration
is aimed to be used in computationally constrained environments, only 20
cepstral coefficients are kept at each window. Although 40 or more cepstral
coefficients are usually used in research papers like [16, 15], it is widely
believed that 15 or more cepstral coefficients will preserve the essential
information in audio signals for classification applications. The resulting
single channel 2D MFCC $\in\mathbb{R}^{20\times 51}$ is finally normalized to
lie between -1 and 1.

###  IV-C Network Structure

For comparison of the proposed layer with the baseline convolutional layer in
terms of accuracy of final speech commands classification, we need to use the
layers in a deep structure which is end-to-end trainable. Since we aim to use
the resulting network in computationally constrained environments, like
embedded devices in robotic applications, we use a lightweight basic model for
the evaluations. Indeed, the model could not dominate the complex state-of-
the-art networks, but it is quite enough for a proof of concept, i.e., showing
superiority of the proposed layer over the baseline convolutional layer. So,
we choose a structure which is similar to LeNet-1 [35]. The original LeNet-1
has two convolutional and one fully connected layers. The network employed
here has the same structure, but the convolutional layers are replaced with
the proposed quadratic SelfONN layer. The codes will also be made freely
available on-line upon publication of the paper.

After each convolutional layer a max pooling layer and a tangent hyperbolic
activation are used. Dropout is also used to prevent over-fitting. The
proposed layer and the resulting network are implemented in PyTorch v1.4.0.

###  IV-D Hyper-parameters

For training of the network, the accuracy over validation set is computed at
each epoch and the training proceeds if both following conditions are met: the
epoch number is less than the maximum allowable epochs (100 in our
simulations) and the accuracy over validation set is not less than the maximum
achieved validation accuracy for the last 10 epochs. The mini-batch size is
set to 50. Based on cross validation on the training set using the CNN
comprised of the baseline convolutional layers, the kernel size is set to
$3\times 3$, and dilation and stride to 1. Padding of size 2 is also applied
to the all sides of the input of the convolutional layers. The number of
neurons (or channels) in two convolutional layers and one fully connected
layer equals to 20, 20, and 10, respectively. Stochastic gradient descent with
momentum = 0.9 and learning rate = 0.01 is used for network training.

###  IV-E Experimental Results

Figures 2 and 3 show the comparison of accuracy versus maximum Taylor series
power ($Q$) over the test set between convolutional, SelfONN, and proposed
layers for GSC and SSC datasets, respectively. As shown, the network comprised
of the proposed layer outperforms the network with convolutional or SelfONN
layers in terms of test accuracy over both datasets. For SelfONN, the values
of $Q=3$ and $Q=4$ result in the best accuracy for GSC and SSC datasets,
respectively. The accuracy of the proposed quadratic SelfONN layer first
increases with increase in $Q$, but the accuracy is saturated when $Q$ further
increases.

As can be inferred from Figures 2 and 3, which is also summarized in Table I,
the accuracy gains for the proposed layer over the baseline convolutional
layer are nearly 4.6% and 2.2% for GSC and SSC datasets, respectively.
Similarly, the accuracy gains of the proposed layer over SelfONN layer are
nearly 2.4% and 0.8% for GSC and SSC datasets, respectively.

Figure 2: Comparison of the convolutional, SelfONN, and proposed layer in
terms of test accuracy for SCR task over GSC dataset.

Figure 3: Comparison of the convolutional, SelfONN, and proposed layer in terms of test accuracy for SCR task over SSC dataset. TABLE I: Performance comparison of the proposed layer with ordinary convolutional and SelfONN layers over GSC and SSC datasets. Dataset | Ordinary convolution | SelfONN | Proposed
---|---|---|---
GSC | 85.2 | 87.6 | 89.8
SSC | 95.5 | 97.1 | 97.9

Figure 4 evaluates the proposed layer based on the deployment complexity. To
this end, the deployment time for a single 1 second speech command is computed
by averaging the deployment time over the whole test set. The simulations and
learning process are run on a machine with a Nvidia GeForce GTX 1080.

As can be seen, the deployment times of the proposed layer are higher than
those of the convolution layer. However, the deployment time of the proposed
layer is quite close to that of SelfONN.

Figure 4: Comparison of the convolutional, SelfONN, and the proposed layers in
terms of deployment time on a Nvidia GeForce GTX 1080 for SCR task over GSC
dataset.

###  IV-F Comparison to the State-of-the-art

According to leader board section of Kaggle challenge on GSC
dataset444https://kaggle.com/c/tensorflow-speech-recognition-challenge , the
first team has achieved an accuracy around 91%. However, they probably used
much more complex networks or ensemble of classifiers. They also used the
original wave signals with full resolution spectrograms or higher number of
cepstral coefficients. Since, the purpose of our work is the proof of concept
of the proposed layer in computationally limited environments, we just used a
very light-weight network, with small resolution inputs and compared the
accuracy with the baseline convolutional layers. Even with these limitations,
the achieved accuracy is quite high. However, it is expected that employing
the proposed layer in deeper networks will result in higher accuracy values. A
similar discussion also holds for SSC dataset.

##  V Conclusion

An extension of the ordinary convolutional layer is proposed which is highly
efficient especially at computationally constrained environments. It is
especially appealing for robotic applications, where the resource limitations
will not allow to use highly deep structures. The performance gap is promising
with the proposed layer for embedded SCR applications. Future work will
involve using the proposed layer in more complex networks in combination with
network squeezing methods like network pruning, distillation, and
sparsification of the fully connected layers.

## References

  * [1] Shinichiro Okada, Yoshiaki Tanaba, Hideyuki Yamauchi, and Shoichi Sato,  “Single-surgeon thoracoscopic surgery with a voice-controlled robot,”  The Lancet, vol. 351, no. 9111, pp. 1249, 1998.
  * [2] Samir A. Elsagheer Mohamed, Allam Shehata Hassanin, and Mohamed Tahar Ben Othman,  “Educational system for the holy quran and its sciences for blind and handicapped people based on google speech api,”  Journal of Software Engineering and Applications, vol. 7, no. 3, pp. 18–30, 2014.
  * [3] Roman A. Solovyev, Alexandr A. Kalinin, Alexander G. Kustov, Dmitry V. Telpukhov, and Vladimir S. Ruhlov,  “FPGA implementation of convolutional neural networks with fixed-point calculations,”  arXiv preprint arXiv:1808.09945, 2018.
  * [4] K.-F. Lee and H.-W. Hon,  “Speaker-independent phone recognition using hidden Markov models,”  IEEE Transactions on Acoustics, Speech, and Signal Processing, vol. 37, no. 11, pp. 1641–1648, 1989.
  * [5] Urmila Shrawankar and Vilas M. Thakare,  “Techniques for feature extraction in speech recognition system: A comparative study,”  arXiv preprint arXiv:1305.1145, 2013.
  * [6] Dan Jurafsky and James H. Martin,  Speech and Language Processing. Vol. 3,  Pearson London London, 2014.
  * [7] Herve A. Bourlard and Nelson Morgan,  Connectionist Speech Recognition: A Hybrid Approach, vol. 247,  Springer Science & Business Media, 2012.
  * [8] James L. McClelland and Jeffrey L. Elman,  “The TRACE model of speech perception,”  Cognitive psychology, vol. 18, no. 1, pp. 1–86, 1986.
  * [9] Geoffrey Hinton, Li Deng, Dong Yu, George E. Dahl, Abdel-rahman Mohamed, Navdeep Jaitly, Andrew Senior, Vincent Vanhoucke, Patrick Nguyen, and Tara N. Sainath,  “Deep neural networks for acoustic modeling in speech recognition: The shared views of four research groups,”  IEEE Signal processing magazine, vol. 29, no. 6, pp. 82–97, 2012\.
  * [10] Yann LeCun, Yoshua Bengio, and Geoffrey Hinton,  “Deep learning,”  nature, vol. 521, no. 7553, pp. 436–444, 2015.
  * [11] Ian Goodfellow, Yoshua Bengio, and Aaron Courville,  Deep Learning,  MIT press, 2016.
  * [12] Felix A. Gers, Jürgen Schmidhuber, and Fred Cummins,  “Learning to forget: Continual prediction with LSTM,”  in 9th International Conference on Artificial Neural Networks. 1999, pp. 850–855, IET.
  * [13] Klaus Greff, Rupesh K. Srivastava, Jan Koutník, Bas R. Steunebrink, and Jürgen Schmidhuber,  “LSTM: A search space odyssey,”  IEEE transactions on neural networks and learning systems, vol. 28, no. 10, pp. 2222–2232, 2016.
  * [14] Guoguo Chen, Carolina Parada, and Georg Heigold,  “Small-footprint keyword spotting using deep neural networks,”  in 2014 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). 2014, pp. 4087–4091, IEEE.
  * [15] Yundong Zhang, Naveen Suda, Liangzhen Lai, and Vikas Chandra,  “Hello edge: Keyword spotting on microcontrollers,”  arXiv preprint arXiv:1711.07128, 2017.
  * [16] Raphael Tang and Jimmy Lin,  “Deep residual learning for small-footprint keyword spotting,”  in 2018 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), 2018, pp. 5484–5488.
  * [17] Seungwoo Choi, Seokjun Seo, Beomjun Shin, Hyeongmin Byun, Martin Kersner, Beomsu Kim, Dongyoung Kim, and Sungjoo Ha,  “Temporal convolution for real-time keyword spotting on mobile devices,”  arXiv preprint arXiv:1904.03814, 2019.
  * [18] Mirco Ravanelli and Yoshua Bengio,  “Speaker recognition from raw waveform with sincnet,”  in IEEE Spoken Language Technology Workshop (SLT), 2018, pp. 1021–1028.
  * [19] Mirco Ravaneli and Yoshua Bengio,  “Interpretable convolutional filters with sincnet,”  arXiv preprint arXiv:1811.09725, 2018.
  * [20] Neil Zeghidour, Nicolas Usunier, Gabriel Synnaeve, Ronan Collobert, and Emmanuel Dupoux,  “End-to-end speech recognition from the raw waveform,”  arXiv preprint arXiv:1806.07098, 2018.
  * [21] Justice Amoh and Kofi M. Odame,  “An Optimized Recurrent Unit for Ultra-Low-Power Keyword Spotting,”  in Proceedings of the ACM on Interactive, Mobile, Wearable and Ubiquitous Technologies. June 2019, vol. 3 of 2, pp. 1–17, Association for Computing Machinery.
  * [22] Théodore Bluche, Maël Primet, and Thibault Gisselbrecht,  “Small-Footprint Open-Vocabulary Keyword Spotting with Quantized LSTM Networks,”  arXiv:2002.10851, Feb. 2020.
  * [23] Serkan Kiranyaz, Turker Ince, Alexandros Iosifidis, and Moncef Gabbouj,  “Operational Neural Networks,”  arXiv:1902.11106, Oct. 2019.
  * [24] Serkan Kiranyaz, Junaid Malik, Habib Ben Abdallah, Turker Ince, Alexandros Iosifidis, and Moncef Gabbouj,  “Self-Organized Operational Neural Networks with Generative Neurons,”  arXiv preprint arXiv:2004.11778, 2020.
  * [25] Serkan Kiranyaz, Turker Ince, Alexandros Iosifidis, and Moncef Gabbouj,  “Progressive operational perceptrons,”  Neurocomputing, vol. 224, pp. 142–154, 2017.
  * [26] Dat Thanh Tran, Serkan Kiranyaz, Moncef Gabbouj, and Alexandros Iosifidis,  “PyGOP: A Python library for Generalized Operational Perceptron algorithms,”  Knowledge-Based Systems, vol. 182, pp. 104801, 2019.
  * [27] DatThanh Tran, Serkan Kiranyaz, Moncef Gabbouj, and Alexandros Iosifidis,  “Progressive Operational Perceptrons with Memory,”  Neurocomputing, vol. 379, pp. 172–181, 2020.
  * [28] Dat Thanh Tran, Serkan Kiranyaz, Moncef Gabbouj, and Alexandros Iosifidis,  “Heterogeneous multilayer generalized operational perceptron,”  IEEE transactions on neural networks and learning systems, vol. 31, no. 3, pp. 710–724, 2019.
  * [29] Dat Thanh Tran, Alexandros Iosifidis, and Moncef Gabbouj,  “Improving efficiency in convolutional neural networks with multilinear filters,”  Neural Networks, vol. 105, pp. 328–339, 2018.
  * [30] François Chollet,  “Xception: Deep learning with depthwise separable convolutions,”  in Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, 2017, pp. 1251–1258.
  * [31] Andrew G. Howard, Menglong Zhu, Bo Chen, Dmitry Kalenichenko, Weijun Wang, Tobias Weyand, Marco Andreetto, and Hartwig Adam,  “Mobilenets: Efficient convolutional neural networks for mobile vision applications,”  arXiv preprint arXiv:1704.04861, 2017.
  * [32] Lukasz Kaiser, Aidan N. Gomez, and Francois Chollet,  “Depthwise separable convolutions for neural machine translation,”  arXiv preprint arXiv:1706.03059, 2017.
  * [33] Georgios Zoumpourlis, Alexandros Doumanoglou, Nicholas Vretos, and Petros Daras,  “Non-linear convolution filters for CNN-based learning,”  in Proceedings of the IEEE International Conference on Computer Vision, 2017, pp. 4761–4769.
  * [34] Pete Warden,  “Speech commands: A dataset for limited-vocabulary speech recognition,”  arXiv preprint arXiv:1804.03209, 2018.
  * [35] Yann LeCun, Léon Bottou, Yoshua Bengio, and Patrick Haffner,  “Gradient-based learning applied to document recognition,”  Proceedings of the IEEE, vol. 86, no. 11, pp. 2278–2324, 1998.


# Mode Delocalization in Disordered Photonic Chern Insulator

Udvas Chattopadhyay  Division of Physics and Applied Physics, School of
Physical and Mathematical Sciences,
Nanyang Technological University, Singapore 637371, Singapore  Sunil Mittal
Department of Electrical and Computer Engineering, The University of Maryland
at College Park, College Park, MD 20742, USA  Mohammad Hafezi  Department of
Electrical and Computer Engineering, The University of Maryland at College
Park, College Park, MD 20742, USA  Y. D. Chong<EMAIL_ADDRESS>Division of
Physics and Applied Physics, School of Physical and Mathematical Sciences,
Nanyang Technological University, Singapore 637371, Singapore  Centre for
Disruptive Photonic Technologies, Nanyang Technological University, Singapore
637371, Singapore

###### Abstract

In disordered two dimensional Chern insulators, a single bulk extended mode is
predicted to exist per band, up to a critical disorder strength; all the other
bulk modes are localized. This behavior contrasts strongly with topologically
trivial two-dimensional phases, whose modes all become localized in the
presence of disorder. Using a tight-binding model of a realistic photonic
Chern insulator, we show that delocalized bulk eigenstates can be observed in
an experimentally realistic setting. This requires the selective use of
resonator losses to suppress topological edge states, and acquiring
sufficiently large ensemble sizes using variable resonator detunings.

##  I Introduction

One of the most intriguing features of the Integer Quantum Hall Effect von
Klitzing et al. (1980); Prange and Girvin (1990); Wysokinski (2000) is the
extraordinary accuracy of quantization in the Hall resistivity of about 1 part
in $10^{9}$. Disorder plays an important role in this phenomenon; without
disorder, the Integer Quantum Hall Effect’s celebrated conductance
quantization plateaus could not exist Prange and Girvin (1990). As the seminal
work of Anderson and co-workers has shown, the effects of disorder are
strongly dependent on the spatial dimensionality Anderson (1958); Abrahams et
al. (1979); Evers and Mirlin (2008). In one dimension, arbitrarily weak
disorder localizes all states, whereas three dimensional systems host
localized states at low energies and extended states at high energies,
separated by a mobility edge. In two dimensions (2D), the effects of disorder
depend on the time-reversal and spin symmetries of the system García-Martínez
et al. (2013). For normal 2D materials, which are in the orthogonal symmetry
class, all states are localized by disorder, similar to the one dimensional
case MacKinnon and Kramer (1981, 1983). For the unitary class, which includes
Integer Quantum Hall systems and other Chern insulators, localization occurs
via a mechanism called “levitation and annihilation”: in the limit of infinite
system size, the introduction of disorder causes all states to localize except
for one state per (topologically nontrivial) band, which remains extended;
with increasing disorder strength, extended states in adjacent bands can move
towards each other and annihilate, producing a transition to a purely-
localized phase Khmelnitskii (1984); Laughlin (1984); Evangelou (1995); Asada
et al. (2002); Onoda et al. (2007).

The inter-plateau longitudinal conductance peaks in the Integer Quantum Hall
Effect constitute the principal experimental evidence for the special bulk
extended states in 2D unitary disordered systems Mares̆ et al. (1999);
Coleridge (2005). In Chern insulators without Landau levels, there is thus far
little experimental evidence for these states, nor for theorized behaviors
such as levitation and annihilation, though many numerical studies have been
performed Evangelou (1995); Asada et al. (2002); Onoda et al. (2007); Xu et
al. (2012); Castro et al. (2015); Qiao et al. (2016). In condensed matter
settings, such experiments are very challenging due to the need to fabricate
large samples with controlled amounts of disorder.

This paper investigates the possibility of using photonics to probe the
localization behavior of 2D Chern insulators. Over the last decade, photonics
has emerged as a versatile setting for realizing phenomena associated with
band topology Ozawa et al. (2019); Kim et al. (2020), including reflectionless
edge transport Raghu and Haldane (2008); Wang et al. (2009); Rechtsman et al.
(2013); Hafezi et al. (2011, 2013); Leykam et al. (2018); Mittal et al.
(2019), topological pumping Kraus et al. (2012); Hu et al. (2015, 2017);
Mittal et al. (2016); Zilberberg et al. (2018), spin and valley Hall edge
statesWu and Hu (2015); Ma and Shvets (2016), Fermi arcs Lu et al. (2013), and
more. Topological phenomena also hold promise for novel device applications in
photonics, such as highly-robust waveguides and delay lines Wang et al.
(2009); Hafezi et al. (2011, 2013), amplifiers Peano et al. (2016), isolators
Zhou et al. (2017), and lasers Harari et al. (2018); Bandres et al. (2018);
Zeng et al. (2020). There are several reasons to consider using photonic
topological insulators to study the localization properties of topological
phases. First, different disorder configurations can be implemented on a
single device by means of optical, thermal, acoustic, or electrical pumps
Villeneuve et al. (1996); Mittal et al. (2016); Leykam et al. (2018), which
should simplify the acquisition of ensembles with many independent disorder
realizations. Second, it is possible to excite any frequency in the band or
band gap via a number of available launching schemes. Third, field
distributions can be observed by near-field imaging or other techniques,
allowing for the accurate and direct determination of quantities such as
localization lengths Hafezi et al. (2013); Mittal et al. (2019). Fourth,
losses can be controllably incorporated into photonic structures El-Ganainy et
al. (2019), which, as we shall see, is helpful for distinguishing the
experimental signatures of bulk delocalization from the effects of topological
edge states. Although photonics has already been extensively employed for the
experimental study of ordinary Anderson localization Segev et al. (2013), it
has never been used to investigate the peculiar localization properties of
bulk states in Chern insulators.

The downside, however, is that photonic Chern insulators must be formed from
deliberately structured photonic media, such as photonic crystals,
metamaterials, coupled resonators, or waveguide arrays Ozawa et al. (2019);
Kim et al. (2020). In some of these platforms, fabrication technologies are
unable to create lattices that are sufficiently large, relative to the unit
cell, for localization studies. Moreover, several platforms exhibit rather
high radiative and material losses; while a loss-induced decay length of (say)
several unit cells might be acceptable for the purposes of demonstrating
topological edge transport, it could complicate localization studies through
the introduction of exponentially decreasing intensity profiles to states that
are supposed to be delocalized.

We focus on a type of 2D photonic Chern insulator consisting of a lattice of
on-chip coupled ring resonators, which has recently been proposed and
implemented Leykam et al. (2018); Mittal et al. (2019). This system is
amenable to theoretical analysis since it can be accurately described by a
tight-binding model Leykam et al. (2018). It has been experimentally realized
using silicon photonics, featuring large lattice sizes of up to $15\times 7$
unit cells and sufficiently low levels of loss that edge transport was
observed over distances of over a dozen unit cells Mittal et al. (2019). The
lattice parameters can be altered by methods such as optical pumping Leykam et
al. (2018); Mittal et al. (2019), so disorder ensembles can be readily
generated in each photonic lattice through spatially inhomogenous pumping,
eliminating the need to fabricate many different samples.

Figure 1: (a) Schematic of a finite lattice of coupled optical resonators,
composed of $N_{x}$ unit cells along $x$ (length) and $N_{y}$ unit cells along
$y$ (width). A loss term $i\gamma$ is added to the sites along the top and
bottom edges. Light is injected into the sites on the first column via
uniformly excited coupling waveguides (orange arrows); there are no coupling
waveguides on the rightmost column, which serves as a closed edge. Inset:
close-up view of one unit cell, showing the selected direction of circulation
within the resonators. (b)–(c) Calculated band structures of the semi-infinite
lattice (infinite length and finite width, with losses omitted), for (b) the
topological phase ($M=0$) and (c) the trivial phase ($M=4J$). Edge states are
shown in blue. In (b), the right-moving (left-moving) edge state is localized
on the bottom (top) edge. In (c), there are edge states localized on the
bottom edge, but these do not span the gap.

Using tight-binding simulations, we show that this platform could be used to
access the delocalization of bulk Chern insulator states and the levitation
and annihilation phenomenon. A clear experimental signature can be achieved
with a lattice size of about $50\times 12$, which is a modest increase
relative to existing experiments Mittal et al. (2019), and far smaller than
the lattices in previous numerical localization studies (which typically
feature sample lengths of $~{}10^{3}$ or more) MacKinnon and Kramer (1981,
1983); Evangelou (1995); Asada et al. (2002); Wang et al. (2015). Silicon-on-
insulator, which is typically used in experiments involving coupled ring
resonators Mittal et al. (2019), has a loss level of $\sim 1\,\textrm{dB/cm}$.
While silicon is the best choice for compactness, it is not particularly a
low-loss platform. Silicon nitride can achieve loss levels of $\sim
1\,\textrm{dB/m}$ Bauters et al. (2011) and we show that such a level of loss
does not affect the key results. The topological edge states of the Chern
insulator tend to conflict with the experimental signature of the delocalized
bulk states, but we find that the former can be suppressed simply by adding
losses to the resonators along the lattice edges. Hence, the photonic lattice
can provide a way to explore the localization behavior of bulk states in
disordered topological insulators, which have thus far resisted in-depth
experimental investigation.

##  II Model

We consider a photonic Chern insulator of a type that has recently been
proposed Leykam et al. (2018) and implemented using silicon photonics Mittal
et al. (2019). As shown in Fig. 1(a), the system consists of an bipartite
square lattice of resonant “site rings” coupled to off-resonant “link rings”.
The site rings occupy two sublattices, denoted by $A$ and $B$, and the link
rings introduce nearest-neighbor and next-nearest-neighbor couplings between
them. Light propagation within the lattice can be decomposed into two
pseudospin sectors (corresponding to clockwise or counterclockwise circulation
in the site rings), which do not interact due to local momentum conservation
at the inter-ring coupling regions Hafezi et al. (2011, 2013). Within each
sector, time reversal symmetry is effectively broken (however, the physical
structure is time reversal symmetric and can be fabricated from ordinary
dielectric materials).

In the absence of disorder and losses, the system is described by the
following tight-binding Hamiltonian Leykam et al. (2018):

$H=(\epsilon+M)\,n_{A}+(\epsilon-M)\,n_{B}+V_{\textrm{nn}}+V_{\textrm{nnn}},$
(1)

where

$\displaystyle\epsilon$ $\displaystyle=2J\cot(\phi/2)$ (2) $\displaystyle
n_{A}$
$\displaystyle=\sum_{r}a_{r}^{\dagger}\,a_{r},\;\;\;n_{B}=\sum_{r}b_{r}^{\dagger}\,b_{r}$
(3) $\displaystyle V_{\textrm{nn}}$
$\displaystyle=W_{1}\sum_{r}\left[a_{r}^{\dagger}(b_{r}+b_{r+x+y})+b_{r}^{\dagger}(a_{r-x}+a_{r-y})\right]$
(4) $\displaystyle V_{\textrm{nnn}}$
$\displaystyle=W_{2}\sum_{r}\Big{[}a_{r}^{\dagger}\sum_{\pm}\hat{a}_{r\pm
y}+b_{r}^{\dagger}\sum_{\pm}b_{r\pm x}\Big{]}$ (5) $\displaystyle W_{1}$
$\displaystyle=J\exp(i\phi/4)\csc(\phi/2),\;\;\;W_{2}=J\csc(\phi/2).$ (6)

Here, $a_{r}$ and $b_{r}$ are the annihilation operators for the $A$ and $B$
sublattices on the unit cell at position $r$ (with $r+x$ denoting the position
one unit cell to the right, etc.), $M$ is a sublattice-dependent resonator
detuning for the site rings, and $J$ and $\phi$ parameterize the couplings
mediated by the link rings. For $\phi\neq 0$, time reversal symmetry is
effectively broken. Each eigenvalue of $H$, denoted by $\delta\nu$,
corresponds to the detuning of a photonic eigenmode relative to a reference
frequency.

Figure 2: Normalized localization length $\xi_{N}$ versus relative detuning
$\delta\nu$ for samples with length $N_{x}=50$, different widths $N_{y}$,
different disorder strengths $W$, and open boundary condition along $y$. A
loss term $-i\gamma$ with $\gamma=J$ is added to the edge sites to suppress
the edge state. The ensemble size used is $500$. Black arrows indicate the
critical behavior ($\xi_{N}$ increasing or constant with $N_{y}$), as shown in
the inset for $W=6$. For large $W$, levitation and annihilation of the
critical states is observed, and the bulk states are completely localized for
$W\geq 10J$. The band gap of the ordered bulk system is highlighted in yellow.

For a given nonzero $\phi$, which depends on the detuning of the link rings
relative to the site rings, the lattice supports both a topological band
insulator phase (a Chern insulator) and a topologically trivial phase (a
normal insulator), depending on the value of $M/J$, the relative detuning
between the site ring sublattices Leykam et al. (2018). In the following, we
take $\phi=\pi$ (i.e., link rings exactly anti-resonant with the site rings),
with $M=0$ for the Chern insulator and $M=4J$ for the normal insulator. In
Fig. 1(b)–(c), their bandstructures are plotted for a disorder-free quasi-one-
dimensional geometry (i.e., a strip that is infinite in $x$, and finite in $y$
with open boundary conditions). For the Chern insulator ($M=0$), the gap is
spanned by chiral edge states that are localized to opposite edges, consistent
with the Chern numbers of $\pm 1$ possessed by the two bulk bands Ozawa et al.
(2019); Kim et al. (2020); Leykam et al. (2018); Mittal et al. (2019). For the
normal insulator ($M=4J$), there are no edge states spanning the gap.

We consider rectangular lattices of length $N_{x}$ and width $N_{y}$. As
indicated in Fig. 1(a), light is injected uniformly into the lattice via
waveguides coupled to the site rings on the left edge (i.e., all site rings
along that edge are excited with equal intensity and phase). The light
returning from the lattice bulk is assumed to outcouple through the same
coupling waveguides.

In an actual experiment, the intensities on individual resonators can be
determined from direct measurements of weak light scattering Hafezi et al.
(2013); Mittal et al. (2019). The site intensities in column $n$ can be
calculated from the frequency-domain Green’s function (which we will obtain
using the method of Kramer and McKinnon MacKinnon and Kramer (1983), which is
based on the recursive Green’s function technique Thouless and Kirkpatrick
(1981); Lewenkopf and Mucciolo (2013); Zhang et al. (2017)). However, the
calculation has a subtle dependence on whether the excitations on column $1$
are mutually coherent or incoherent. First, consider the incoherent case, in
which different sites bear no fixed phase relationship with one another. The
total intensity on column $n$ is
$\mathrm{Tr}\left(\mathcal{G}^{\dagger}\mathcal{G}\right)$, where
$\mathcal{G}$ is the Green’s function matrix beween sites on column $1$ and
sites on column $n$. We can thus define

$T=\frac{1}{N_{y}}\,\mathrm{Tr}\left(\mathcal{G}^{\dagger}\mathcal{G}\right),$
(7)

which is the standard definition of the transmittance as used in the
electronic transport literature for determining the conductance of a sample
MacKinnon and Kramer (1983). On the other hand, in photonic experiments the
inputs typically originate from a single laser source with strong spatial
coherence. In that case, the effective transmittance in column $n$ is

$T_{c}=\frac{\langle\psi_{1}|\mathcal{G}^{\dagger}\mathcal{G}|\psi_{1}\rangle}{\langle\psi_{1}|\psi_{1}\rangle},$
(8)

where $|\psi_{1}\rangle$ is the input vector (e.g., $[1,1,\dots,1]$ if all
input waveguides have the same phase).

Anderson-type disorder is introduced into the lattice in the form of a random
detuning on each site ring, drawn independently from a uniform distribution
over $[-W/2,W/2]$, where $W$ is a tunable disorder strength parameter. In an
actual photonic lattice, such disorder is introduced in part by unavoidable
fabrication imperfections; a previous experiment found this to be on the order
of the coupling strength $J$ (the basic “energy” scale of the tight-binding
model) Mittal et al. (2019). Additional disorder can be introduced through
spatially inhomogenous pumping Villeneuve et al. (1996); Mittal et al. (2016);
Leykam et al. (2018). We assume that an ensemble of independent disorder
configurations can thus be achieved.

We can then compute the disorder average $\langle\log(T)\rangle$ for the
incoherent input case. This is related to the localization length $\xi$ by
Slevin et al. (2001); Paulin and Carpentier (2012)

$\langle\log(T)\rangle\propto-\frac{2n}{\xi}.$ (9)

Therefore, $\xi$ can be extracted from a linear fit of $\langle\log(T)\rangle$
against $n$. For coherent inputs, we substitute $T$ with $T_{c}$ in Eq. (9).
We find numerically that although the values of $\langle\log T\rangle$ and
$\langle\log T_{c}\rangle$ are generally different, both cases yield the same
fitted value of $\xi$, as shown in Appendix A. Hence, the use of spatially
coherent inputs in a photonic experiment is consistent with the standard
definition of the localization length, based on Eq. (7), which is what our
subsequent numerical results are based on.

We then define the normalized localization length MacKinnon and Kramer (1981)

$\xi_{N}(\delta\nu,W)=\frac{\xi}{N_{y}},$ (10)

which depends implicitly on the operating frequency detuning $\delta\nu$, as
well as the disorder strength $W$. In the localized regime, $\xi_{N}$
decreases with $N_{y}$ and vanishes in the $N_{y}\to\infty$ limit. When
extended states are present, $\xi_{N}$ increases with $N_{y}$ and diverges in
the $N_{y}\to\infty$ limit. For the critical states at a mobility edge,
$\xi_{N}$ approaches a finite constant in as $N_{y}\to\infty$. Hence, for
given $\delta\nu$ and $W$, we can detect the presence of extended states by
finding how $\xi_{N}(\delta\nu)$ varies with $N_{y}$ Xu et al. (2012);
Evangelou (1995); Su et al. (2016); Wang et al. (2015); Castro et al. (2016).

In numerical studies of bulk localization in 2D lattices, it has been
conventional to take periodic boundary conditions along the upper and lower
edges, so that the sample forms an edgeless waveguide MacKinnon and Kramer
(1981, 1983); Asada et al. (2002); Wang et al. (2015). The reason for doing
this is that Chern insulators and other 2D topological insulators host
topological edge states that are robust against disorder and extended along
the edge, which can interfere with the signature of bulk localization of
delocalization. However, since we aim to investigate the feasibility of a
realistic on-chip photonic experiment, it is not appropriate for us to impose
periodic boundary conditions. Instead, we introduce losses to the sites
residing at the edge of the system, as shown by the gray outlines in Fig.
1(a). The losses are modelled as an imaginary contribution to the detuning,
$-i\gamma$, where we take $\gamma=J$ (later, we will also study the effects of
smaller losses on all the other sites). In experiments, losses can be
deliberately introduced by using adding lossy materials or claddings or
additional scattering defects to the resonators El-Ganainy et al. (2019); Gao
et al. (2016). It would be desirable to avoid significantly altering the real
detuning of the edge resonators while doing so; otherwise, the edge states may
simply be shifted onto adjacent rows further into the bulk.

##  III Results

Fig. 2 shows the normalized localization length $\xi_{N}$ versus the source
frequency detuning $\delta\nu$, for different disorder strengths $W$ and
different lattice widths $N_{y}$. The lattices are in the Chern insulator
phase ($M=0$), and have fixed length $N_{x}=50$. At the upper and lower edges,
we impose open (Dirichlet-like) boundary conditions. As mentioned, losses are
added to the sites on the upper and lower edges to suppress the topological
edge state; the losses on the other resonators are assumed to be negligible.

For most values of $\delta\nu$, we find that $\xi_{N}$ decreases with
increasing $N_{y}$. As explained in the previous section, this indicates that
the eigenmodes at these frequencies are localized by the disorder. However, in
Fig. 2(a)–(d), representing disorder strengths up to $W=8J$, there are two
regions on each side of the band gap where $\xi_{N}$ is constant or increases
with $N_{y}$. This is evidence for the bulk extended states believed to exist
in disordered Chern insulators Evangelou (1995); Asada et al. (2002); Onoda et
al. (2007); Xu et al. (2012); Castro et al. (2015); Qiao et al. (2016). In
agreement with theoretical predictions, there is only a narrow range of
frequencies within each band where this phenomenon occurs.

The signature of delocalization is more easily observable if the ensemble size
is large. In Fig. 2, each data point is averaged from an ensemble of $500$
independent disorder realizations. It would be experimentally unfeasible to
fabricate one physical sample for each disorder realization, but it should be
possible to implement disorder dynamically via spatially inhomogenous and
actively-switchable optical, thermal, acoustic, or electrical pumps Villeneuve
et al. (1996); Mittal et al. (2016); Leykam et al. (2018). This would allow
numerous independent disorder realizations to be generated with a few physical
samples (corresponding to different values of $N_{y}$). Furthermore, according
to our simulations, delocalization may still be observable for ensemble sizes
of as low as $50$.

Figure 3: (a) Chern number of the lower band versus disorder strength $W$ for
a $32\times 32$ lattice (averaged over $200$ disorder realizations). (b)
Critical energy $\delta\nu_{c}$ (where the extended state is located) versus
$W$ for the realistic sample corresponding to Fig. 2. Each shaded region
indicates the detuning range over which $\xi_{N}$ does not decrease with
$N_{y}$, and the dots and solid lines indicate the mid-point of the range. (c)
Normalized localization length versus $\delta\nu$ for the hermitian quasi-1D
system ($N_{x}=10000$) with periodic boundary along $y$ at a disorder strength
$W=8J$. (d) Critical energy $\delta\nu_{c}$ versus $W$ for semi infinite
lattice corresponding to (c). Anderson localization transition by annihilation
can be observed at $W=9J$.

Fig. 3(a) shows the Chern number calculated for the lower band versus the
disorder strength $W$ Zhang et al. (2013); Castro et al. (2015). The Chern
number is computed by summing the Berry flux through each plaquette of the
discretized Brillouin zone following Fukui’s method Fukui et al. (2005), using
a computationally efficient coupling-matrix method to calculate the Berry flux
Zhang et al. (2013); Castro et al. (2015). With increasing $W$, the Chern
number is quantized to unity until $W\approx 8J$, after which it decreases to
zero. This is consistent with the results shown in Fig. 2, where the signature
of the delocalized bulk mode, a property of the Chern insulator, disappears
between $W=8J$ and $W=10J$. In Fig. 3(b), we plot the range of critical
frequencies (i.e., the frequency range over which $\xi_{N}$ does not decrease
with $N_{y}$) versus $W$. As the system transitions from Chern insulator to
normal insulator, the critical frequencies shift toward one another, and for
$W\gtrsim 9J$ become impossible to distinguish over the statistical noise.
Hence, a photonic lattice can provide evidence for the long-standing
theoretical prediction that topologically insulating behavior is destroyed by
strong disorder through the levitation and annihilation of the delocalized
bulk states Khmelnitskii (1984); Laughlin (1984); Evangelou (1995); Asada et
al. (2002); Onoda et al. (2007).

In Fig. 3(c), we present simulation results for an idealized and
experimentally unrealistic lattice with periodic boundary conditions on the
upper and lower edges, and a much greater length of $N_{x}=10^{4}$.
Delocalization is observed at around the same frequencies as in Fig. 2(d).
This shows that the results in Fig. 2, which were obtained for experimentally
realistic open boundary conditions and much shorter length $N_{x}=50$,
accurately capture the frequencies at which the delocalized bulk states are
supposed to occur. In Fig. 3(d), we plot the range of critical frequencies
versus disorder strength $W$ using the idealized (periodic boundary conditions
and large $N_{x}$) lattice. It displays the same levitation and annihilation
behavior as in Fig. 3(b). The large error bar near the annihilation point
($W\approx 9J$) is due to the fact that the two humps seen in Fig. 3(c) merge
into a single large hump near the annihilation point, for the range of $N_{y}$
considered here. In Fig. 3(b), this was not observed due to the smaller
lattice size and the effects of the open boundary conditions.

Figure 4: Plots of normalized localization length versus relative detuning
$\delta\nu$ for different scenarios. (a) Chern insulator ($M=0$) with disorder
strength $W=6J$, with the losses on the upper and lower edge sites omitted.
Transmission now peaks in the bulk gap (yellow region) due to the topological
edge states, and the signature of bulk delocalization cannot be discerned. (b)
Chern insulator ($M=0$) with disorder strength $W=6J$, and with small losses
$-i\gamma^{\prime}$ on the non-edge sites. Here we take
$\gamma^{\prime}=10^{-3}J$. As before, the sites on the upper and lower edges
have loss $-iJ$, and all other parameters are the same as in Fig. 2. Bulk
delocalization can be observed at $\delta\nu\approx-2J$. For clarity, only
negative values of $\delta\nu$ are plotted. (c)–(d) Normal insulator ($M=4J$)
with disorder strengths of (c) $4J$ and (d) $6J$ and length $N_{x}=50$.

To verify that the loss on the edge sites is necessary for the observation of
bulk delocalization, in Fig. 4(a) we plot $\xi_{N}$ versus $\delta\nu$ with
these losses omitted. In this case, high transmission is observed in the
frequency range corresponding to the bulk gap, due to transport by the now-
unsuppressed topological edge states. The frequency range over which $\xi_{N}$
is constant or increasing with $N_{y}$ appears to occur near the center of the
bulk gap, rather than within each band as in Fig. 2.

In Fig. 4(b), we show that the delocalization signature can still be observed
when there is weak but nonzero losses on the other resonators. In real
photonic structures, some material and radiative loss is always present. Here,
we assign each non-edge site a loss of $-i\gamma^{\prime}$, where
$\gamma^{\prime}=10^{-3}J$. This level of losses is consistent what can be
achieved experimentally. For example, the silicon-on-insulator implementation
of the model in Ref. Mittal et al., 2019 had a loss of $\gamma^{\prime}\approx
0.03J$ and a silicon nitride platform Bauters et al. (2011) can achieve a
loss-level smaller by a factor of around $1/100$. Therefore, a loss level of
$\gamma^{\prime}=10^{-3}J$ should be achievable, and it should hence be
possible to observe a frequency region in which $\xi_{N}$ is approximately
constant with $N_{y}$.

Finally, to verify the topological origin of the bulk extended state, Fig.
4(c)–(d) shows the results for $M=4J$, for which the lattice is in its normal
insulator phase. In this case, there is no clear range of frequencies in which
$\xi_{N}$ is constant or increasing with $N_{y}$.

##  IV Discussion

We have proposed an experimentally feasible way to probe the localization
behavior of disordered Chern insulators using a recently-developed photonic
platform Mittal et al. (2019). Using a realistic tight-binding model Leykam et
al. (2018), we showed that it should be possible to observe the existence of
extended states in the bulk bands, a characteristic feature of disordered
Chern insulators, as well as the levitation and annihilation of these extended
states under increasing disorder Khmelnitskii (1984); Laughlin (1984);
Evangelou (1995); Asada et al. (2002); Onoda et al. (2007). The required
system sizes, disorder strengths, and loss levels are all in the
experimentally accessible range.

There may be other ways to use photonic Chern insulators to probe the
interplay of disorder and topological phases, such as level-spacing statistics
Evers and Mirlin (2008); Castro et al. (2016). Apart from observing
localization lengths from averaged intensity measurements, it may also be
possible to probe the system using the Wigner time delay, which is insensitive
to losses and has previously been used to establish the extended nature of
photonic topological edge states Mittal et al. (2014); however, we have thus
far been unable to find a clear signature of bulk state delocalization in
Wigner time delay statistics. It may also be interesting to explore the
effects of adding optical gain to such a system, to determine whether bulk
extended modes could be observed through their promotion into lasing modes.

Finally, although our proposal has focused on photonic resonator lattices,
similar ideas can be generalized to other photonic lattices, and more broadly
to other bosonic systems, such as acoustic, electrical or mechanical lattices
Yang et al. (2015); Lee et al. (2018); Süsstrunk and Huber (2015).

###### Acknowledgements.

We thank Prof. Caio Lewenkopf for helpful discussions. CYD and UC acknowledge
support from the Singapore Ministry of Education Tier 3 grant MOE2016-T3-1-006
and Tier 1 grant RG187/18. MF and SM acknowledge support from the US Air Force
Office of Scientific Research (AFOSR) Multidisciplinary University Research
Ini-tiative (MURI) grants FA95501610323 and FA95502010223, the US Office of
Naval Research (ONR) MURI grant N00014-20-1-2325, and the National Science
Foundation grant PHY1820938.

Figure 5: (a) Plots of $\langle\log(T)\rangle$ and $\langle\log(T_{c})\rangle$
versus column index $n$ for samples of length $N_{x}=50$, width $N_{y}=10$,
and disorder strength $W=6J$, averaging over 500 disorder realizations. The
straight lines are linear least squares fits. The fitted slopes are almost
identical, so incoherent and coherent excitation give the same localization
length estimate. (b) Normalized localization length $\xi_{N}$ extracted using
$T$ and $T_{c}$, plotted against relative detuning $\delta\nu/J$ for disorder
strength $W=6J$. All other parameters are the same as in Fig. 2.

##  Appendix A Coherent versus Incoherent Excitation

Let $\mathcal{G}_{ba}$ denote the Green’s function matrix element between site
$a$ on the input column $1$ to site $b$ on column $n$. If the input on site
$a$ is $\psi^{(\mathrm{1})}_{a}$, the complex wave amplitude on site $b$ is
$\psi_{b}^{(n)}=\sum_{a}\mathcal{G}_{ba}\psi^{(1)}_{a}$, and the intensity on
that site is

$I_{b}^{(n)}=\left|\psi_{b}^{(n)}\right|^{2}=\sum_{ac}G_{ba}G_{bc}^{*}\psi^{(\mathrm{in})}_{a}\psi^{(\mathrm{in})*}_{c}.$
(11)

For uniform incoherent excitation (i.e., no fixed phase relationship between
different input sites), we take an average over an ensemble of input wave
amplitudes with equal magnitude $\sqrt{I_{\mathrm{in}}/N_{y}}$ and random
phases. Averaging over the ensemble gives the mean intensity

$\left\langle
I_{b}^{(n)}\right\rangle=\frac{I_{\mathrm{in}}}{N_{y}}\sum_{a}\left|G_{ba}\right|^{2},$
(12)

Hence, the total intensity in column $n$ is

$I^{(n)}=\sum_{b}\left\langle
I_{b}^{(n)}\right\rangle=\frac{I_{\mathrm{in}}}{N_{y}}\;\mathrm{Tr}\left(\mathcal{G}_{n}^{\dagger}\mathcal{G}_{n}\right).$
(13)

Normalizing by the total input intensity $I_{\mathrm{in}}$ yields the
transmittance formula given in Eq. (7). On the other hand, if the input is
spatially coherent, we do not perform an average over the random phases, and
the intensity depends explicitly on the input vector, leading to Eq. (8).

As shown in Fig. 5(a), $\langle\log(T)\rangle$ and $\langle\log(T_{c})\rangle$
exhibit the same scaling with distance $n$. Hence, the same estimate for the
localization length $\xi_{N}$ is obtained regardless of whether $T$ or $T_{c}$
is used in Eq. (9). As an example, Fig. 5(b) shows the plot of $\xi_{N}$
versus $\delta\nu$, similar to Fig. 2, for disorder strength $W=6J$.
Estimating the localization length using either $T$ (incoherent excitation) or
$T_{c}$ (coherent excitation) gives essentially the same results.

## References

  * von Klitzing et al. (1980) K. von Klitzing, G. Dorda, and M. Pepper, Phys. Rev. Lett. 45, 494 (1980).
  * Prange and Girvin (1990) R. E. Prange and S. M. Girvin, _The Quantum Hall Effect_ (Springer-Verlag, 1990).
  * Wysokinski (2000) K. I. Wysokinski, European Journal of Physics 21, 535 (2000).
  * Anderson (1958) P. W. Anderson, Phys. Rev. 109, 1492 (1958).
  * Abrahams et al. (1979) E. Abrahams, P. W. Anderson, D. C. Licciardello, and T. V. Ramakrishnan, Phys. Rev. Lett. 42, 673 (1979).
  * Evers and Mirlin (2008) F. Evers and A. D. Mirlin, Rev. Mod. Phys. 80, 1355 (2008).
  * García-Martínez et al. (2013) N. A. García-Martínez, A. G. Grushin, T. Neupert, B. Valenzuela, and E. V. Castro, Phys. Rev. B 88, 245123 (2013).
  * MacKinnon and Kramer (1981) A. MacKinnon and B. Kramer, Phys. Rev. Lett. 47, 1546 (1981).
  * MacKinnon and Kramer (1983) A. MacKinnon and B. Kramer, Zeitschrift für Physik B Condensed Matter 53, 1 (1983).
  * Khmelnitskii (1984) D. E. Khmelnitskii, Phys. Lett. A 106, 182 (1984).
  * Laughlin (1984) R. B. Laughlin, Phys. Rev. Lett. 52, 2304 (1984).
  * Evangelou (1995) S. N. Evangelou, Phys. Rev. Lett. 75, 2550 (1995).
  * Asada et al. (2002) Y. Asada, K. Slevin, and T. Ohtsuki, Phys. Rev. Lett. 89, 256601 (2002).
  * Onoda et al. (2007) M. Onoda, Y. Avishai, and N. Nagaosa, Phys. Rev. Lett. 98, 076802 (2007).
  * Mares̆ et al. (1999) J. J. Mares̆, J. Kris̆tofik, and P. Hubík, Phys. Rev. Lett. 82, 4699 (1999).
  * Coleridge (2005) P. T. Coleridge, in _Encyclopedia of Condensed Matter Physics_ , edited by F. Bassani, G. L. Liedl, and P. Wyder (Elsevier, 2005), p. 248, ISBN 978-0-12-369401-0.
  * Xu et al. (2012) Z. Xu, L. Sheng, D. Y. Xing, E. Prodan, and D. N. Sheng, Phys. Rev. B 85, 075115 (2012).
  * Castro et al. (2015) E. V. Castro, M. P. López-Sancho, and M. A. H. Vozmediano, Phys. Rev. B 92, 085410 (2015).
  * Qiao et al. (2016) Z. Qiao, Y. Han, L. Zhang, K. Wang, X. Deng, H. Jiang, S. A. Yang, J. Wang, and Q. Niu, Phys. Rev. Lett. 117, 056802 (2016).
  * Ozawa et al. (2019) T. Ozawa, H. M. Price, A. Amo, N. Goldman, M. Hafezi, L. Lu, M. C. Rechtsman, D. Schuster, J. Simon, O. Zilberberg, et al., Rev. Mod. Phys. 91, 015006 (2019).
  * Kim et al. (2020) M. Kim, Z. Jacob, and J. Rho, Light: Science & Applications 9, 130 (2020).
  * Raghu and Haldane (2008) S. Raghu and F. D. M. Haldane, Phys. Rev. A 78, 033834 (2008).
  * Wang et al. (2009) Z. Wang, Y. Chong, J. D. Joannopoulos, and M. Soljačić, Nature 461, 772 (2009).
  * Rechtsman et al. (2013) M. C. Rechtsman, J. M. Zeuner, Y. Plotnik, Y. Lumer, D. Podolsky, F. Dreisow, S. Nolte, M. Segev, and A. Szameit, Nature 496, 196 (2013).
  * Hafezi et al. (2011) M. Hafezi, E. A. Demler, M. D. Lukin, and J. M. Taylor, Nat. Phys. 7, 907 (2011).
  * Hafezi et al. (2013) M. Hafezi, S. Mittal, J. Fan, A. Migdall, and J. M. Taylor, Nat. Photon. 7, 1001 (2013).
  * Leykam et al. (2018) D. Leykam, S. Mittal, M. Hafezi, and Y. D. Chong, Phys. Rev. Lett. 121, 023901 (2018).
  * Mittal et al. (2019) S. Mittal, V. V. Orre, D. Leykam, Y. D. Chong, and M. Hafezi, Phys. Rev. Lett. 123, 043201 (2019).
  * Kraus et al. (2012) Y. E. Kraus, Y. Lahini, Z. Ringel, M. Verbin, and O. Zilberberg, Phys. Rev. Lett. 109, 106402 (2012).
  * Hu et al. (2015) W. Hu, J. C. Pillay, K. Wu, M. Pasek, P. P. Shum, and Y. D. Chong, Phys. Rev. X 5, 011012 (2015).
  * Hu et al. (2017) W. Hu, H. Wang, P. P. Shum, and Y. D. Chong, Phys. Rev. B 95, 184306 (2017).
  * Mittal et al. (2016) S. Mittal, S. Ganeshan, J. Fan, A. Vaezi, and M. Hafezi, Nat. Photon. 10, 180 (2016).
  * Zilberberg et al. (2018) O. Zilberberg, S. Huang, J. Guglielmon, M. Wang, K. P. Chen, Y. E. Kraus, and M. C. Rechtsman, Nature 553, 59 (2018).
  * Wu and Hu (2015) L.-H. Wu and X. Hu, Phys. Rev. Lett. 114, 223901 (2015).
  * Ma and Shvets (2016) T. Ma and G. Shvets, New Journal of Physics 18, 025012 (2016).
  * Lu et al. (2013) L. Lu, L. Fu, J. D. Joannopoulos, and M. Soljačić, Nat. Photon. 7, 294 (2013).
  * Peano et al. (2016) V. Peano, M. Houde, F. Marquardt, and A. A. Clerk, Phys. Rev. X 6, 041026 (2016).
  * Zhou et al. (2017) X. Zhou, Y. Wang, D. Leykam, and Y. D. Chong, New J. Phys. 19, 095002 (2017).
  * Harari et al. (2018) G. Harari, M. A. Bandres, Y. Lumer, M. C. Rechtsman, Y. D. Chong, M. Khajavikhan, D. N. Christodoulides, and M. Segev, Science 359, eaar4003 (2018).
  * Bandres et al. (2018) M. A. Bandres, S. Wittek, G. Harari, M. Parto, J. Ren, M. Segev, D. N. Christodoulides, and M. Khajavikhan, Science 359, eaar4005 (2018).
  * Zeng et al. (2020) Y. Zeng, U. Chattopadhyay, B. Zhu, B. Qiang, J. Li, Y. Jin, L. Li, A. G. Davies, E. H. Linfield, B. Zhang, et al., Nature 578, 246 (2020).
  * Villeneuve et al. (1996) P. R. Villeneuve, S. Fan, and J. D. Joannopoulos, Phys. Rev. B 54, 7837 (1996).
  * El-Ganainy et al. (2019) R. El-Ganainy, M. Khajavikhan, D. N. Christodoulides, and S. K. Ozdemir, Comm. Phys. 2, 37 (2019).
  * Segev et al. (2013) M. Segev, Y. Silberberg, and D. N. Christodoulides, Nat. Photon. 7, 197 (2013).
  * Wang et al. (2015) C. Wang, Y. Su, Y. Avishai, Y. Meir, and X. R. Wang, Phys. Rev. Lett. 114, 096803 (2015).
  * Bauters et al. (2011) J. F. Bauters, M. J. R. Heck, D. John, D. Dai, M.-C. Tien, J. S. Barton, A. Leinse, R. G. Heideman, D. J. Blumenthal, and J. E. Bowers, Opt. Express 19, 3163 (2011).
  * Thouless and Kirkpatrick (1981) D. J. Thouless and S. Kirkpatrick, J. Phys. C: Solid State Phys. 14, 235 (1981).
  * Lewenkopf and Mucciolo (2013) C. H. Lewenkopf and E. R. Mucciolo, Journal of Computational Electronics 12, 203 (2013), ISSN 15698025, eprint 1304.3934.
  * Zhang et al. (2017) S.-H. Zhang, W. Yang, and K. Chang, Phys. Rev. B 95, 075421 (2017).
  * Slevin et al. (2001) K. Slevin, P. Markoš, and T. Ohtsuki, Phys. Rev. Lett. 86, 3594 (2001).
  * Paulin and Carpentier (2012) G. Paulin and D. Carpentier, New J. Physics 14, 023026 (2012).
  * Su et al. (2016) Y. Su, C. Wang, Y. Avishai, Y. Meir, and X. R. Wang, Scientific Reports 6, 33304 (2016).
  * Castro et al. (2016) E. V. Castro, R. de Gail, M. P. López-Sancho, and M. A. H. Vozmediano, Phys. Rev. B 93, 245414 (2016).
  * Gao et al. (2016) F. Gao, Z. Gao, X. Shi, Z. Yang, X. Lin, H. Xu, J. D. Joannopoulos, M. Soljačić, H. Chen, L. Lu, et al., Nature Communications 7, 11619 (2016).
  * Zhang et al. (2013) Y.-F. Zhang, Y.-Y. Yang, Y. Ju, L. Sheng, R. Shen, D.-N. Sheng, and D.-Y. Xing, Chinese Physics B 22, 117312 (2013).
  * Fukui et al. (2005) T. Fukui, Y. Hatsugai, and H. Suzuki, Phys. Soc. Jpn. 74, 1674–1677 (2005).
  * Mittal et al. (2014) S. Mittal, J. Fan, S. Faez, A. Migdall, J. M. Taylor, and M. Hafezi, Phys. Rev. Lett. 113, 087403 (2014).
  * Yang et al. (2015) Z. Yang, F. Gao, X. Shi, X. Lin, Z. Gao, Y. Chong, and B. Zhang, Phys. Rev. Lett. 114, 114301 (2015).
  * Lee et al. (2018) C. H. Lee, S. Imhof, C. Berger, F. Bayer, J. Brehm, L. W. Molenkamp, T. Kiessling, and R. Thomale, Communications Physics 1, 39 (2018).
  * Süsstrunk and Huber (2015) R. Süsstrunk and S. D. Huber, Science 349, 47 (2015).