# Algorithm and Hardware Co-Design of Energy-Efficient LSTM Networks for Video
Recognition with Hierarchical Tucker Tensor Decomposition

Yu Gong, Miao Yin, Lingyi Huang, Chunhua Deng, Yang Sui, and Bo Yuan Yu Gong,
Miao Yin, Lingyi Huang, Yang Sui and Bo Yuan are with the Department of
Electrical and Computer Engineering, Rutgers University, Piscataway, NJ,
08854. E-mail<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>Chunhua Deng is currently with ScaleFlux Inc. This work was done when the
author was with Rutgers University. E-mail<EMAIL_ADDRESS>Yang
Sui contributed to this work but was accidentally omitted in the published
version.

###### Abstract

Long short-term memory (LSTM) is a type of powerful deep neural network that
has been widely used in many sequence analysis and modeling applications.
However, the large model size problem of LSTM networks make their practical
deployment still very challenging, especially for the video recognition tasks
that require high-dimensional input data. Aiming to overcome this limitation
and fully unlock the potentials of LSTM models, in this paper we propose to
perform algorithm and hardware co-design towards high-performance energy-
efficient LSTM networks. At algorithm level, we propose to develop fully
decomposed hierarchical Tucker (FDHT) structure-based LSTM, namely FDHT-LSTM,
which enjoys ultra-low model complexity while still achieving high accuracy.
In order to fully reap such attractive algorithmic benefit, we further develop
the corresponding customized hardware architecture to support the efficient
execution of the proposed FDHT-LSTM model. With the delicate design of memory
access scheme, the complicated matrix transformation can be efficiently
supported by the underlying hardware without any access conflict in an on-the-
fly way. Our evaluation results show that both the proposed ultra-compact
FDHT-LSTM models and the corresponding hardware accelerator achieve very high
performance. Compared with the state-of-the-art compressed LSTM models, FDHT-
LSTM enjoys both order-of-magnitude reduction (more than $1000\times$) in
model size and significant accuracy improvement (0.6% to 12.7%) across
different video recognition datasets. Meanwhile, compared with the state-of-
the-art tensor decomposed model-oriented hardware TIE, our proposed FDHT-LSTM
architecture achieve $2.5\times$, $1.46\times$ and $2.41\times$ increase in
throughput, area efficiency and energy efficiency, respectively on LSTM-
Youtube workload. For LSTM-UCF workload, our proposed design also outperforms
TIE with $1.9\times$ higher throughput, $1.83\times$ higher energy efficiency
and comparable area efficiency.

###### Index Terms:

Long Short-term Memory (LSTM), Hardware Architecture, Tensor Decomposition,
Hierarchical Tucker (HT), Video Recognition.

##  1 Introduction

Long short-term memory (LSTM) networks are popularly used in many real-world
sequential analysis and processing tasks. Thanks to their inherent strong
capability of capturing the temporal dependency and modeling the sequential
correlation, the state-of-the-art LSTMs have achieved unprecedented success in
many important temporal sequence-involved applications, such as video
recognition [7], speech recognition [9] and natural language processing (NLP)
[3].

Despite their current widespread adoptions, the large model sizes of LSTM
networks make the practical deployment still very challenging. To be specific,
because the input data of many real-world applications, e.g., video processing
and NLP, naturally exhibit high dimensionality, the corresponding input-to-
hidden weight matrices of LSTMs are typically extremely large. For instance,
as analyzed in [20], even with small-size (256 states) hidden layer, the LSTM
evaluated on UCF11 video recognition dataset [16] already requires more than
50 million weights. Evidently, such ultra-large model size poses a series of
severe challenges for the efficient deployment of LSTMs, including but not
limited to high memory footprint, long processing time, low energy efficiency,
and insufficient accuracy incurred by high difficulty of training.

To address these challenges, various model compression methods, such as
pruning and quantization, have been proposed to build compact LSTMs. Consider
the compression ratio provided by these methods are typically limited and may
be insufficient for the very large LSTMs, tensor decomposition, as a low-rank
approximation approach, have attracted a lot of attention towards high-
performance compression of LSTM models. In general, tensor decomposition aims
to factorize a large-size tensor to a series of small tensor cores with
maintaining small approximation error. After performing tensor decomposition,
a very elegant mathematical phenomenon is that the resulting space and time
complexity can be significantly reduced, sometimes even with order-of-
magnitude reduction. Motivated by this attractive property, recently several
different tensor decomposition approaches, such as tensor train (TT), tensor
ring (TR) and block-term (BT), have been adopted to build compact LSTMs models
[6] [17] [21] for various temporal data analysis tasks, e.g., video
recognition and long-term dynamic forecasting. Compared with the original
uncompressed large-size LSTMs, these tensor decomposed models exhibit much
smaller memory footprints and competitive classification/prediction accuracy.

Although the state-of-the-art tensor decomposed LSTMs demonstrate their
promising potentials, they are still facing two inherent limitations. First,
the existing approaches can only factorize the input-to-hidden layers instead
of the entire model, thereby limiting the overall compression performance.
Notice that, as will be shown in Section 3, the straightforward decomposition
on the hidden-to-hidden layers cannot solve this problem due to the
significant accuracy degradation. Second, the tensor decomposition approaches
adopted in the existing LSTM compression works have inherent limitations. To
be specific, TT decomposition requires that the border tensor cores be rank-1
tensors, thereby directly limiting the overall representation power. Also,
using BT decomposition brings extra flatten and permutation operations, which
causes significant computation overhead for the resulting BT-LSTM models. More
importantly, as will be analyzed in Section 3, the existing adopted tensor
decomposition approaches (TT, TR and BT) in model compression, theoretically,
do not provide the best space complexity reduction. Consequently, they are not
the ideal solutions for building tensor factorized high-accuracy ultra-compact
LSTM models.

To overcome these limitations and develop high-performance LSTM models, in
this paper 111This paper is the extension of the authors’ prior work [22]. we
propose to leverage Hierarchical Tucker (HT) [10] technique, an under-explored
yet powerful tensor tool, to fully decompose the LSTM networks. The resulting
compressed model, namely FDHT-LSTM, enjoys ultra-low complexity with still
maintaining high accuracy. Then, to fully reap the benefits brought by the
proposed efficient compression approach, we further develop the corresponding
hardware architecture to accelerate the execution of FDHT-LSTM models.
Overall, the contributions of this paper are summarized as follows:

  * •

At the algorithm level, we propose to impose fully decomposed hierarchical
Tucker (FDHT) structure on the LSTM networks to achieve ultra-high compression
performance. The proposed FDHT structure contains two main features (as shown
in Figure 1). First, the underlying LSTM is built via using HT decomposition,
a powerful approach that can properly capture and model the correlation and
structure in high-dimensional data. Second, the entire LSTM model, instead of
one or few component layers, is fully decomposed in the HT format in a
homogeneous way. In other words, the low-rank correlation among different
component layers are fully exploited, and thereby bringing order-of-magnitude
reduction in model size while still achieving very high accuracy.

  * •

At the hardware level, we design and implement the corresponding FDHT-LSTM
hardware architecture to fully leverage the algorithmic benefits. In order to
support various complicated matrix transformations introduced by the
computation flow of FDHT-LSTM, we first propose to utilize 2-D SRAM array to
properly map the high-order tensor transpose on the physical 2-D memory
component. Based on this philosophy, we further propose novel write and read
access schemes to the 2-D SRAM array, thereby enabling conflict-free memory
access with high flexibility and reconfigurability to different workloads.

  * •

Evaluation results demonstrate the superior performance of our proposed FDHT-
LSTM and the corresponding hardware accelerator. Experiments show that the
FDHT-LSTM models can use very few parameters (less than 10,000 weights) to
achieve very high accuracy across different video recognition datasets.
Compared with the state-of-the-art compressed LSTM models, FDHT-LSTM enjoys
both order-of-magnitude reduction (more than $1000\times$) in model size and
significant accuracy improvement (0.6% to 12.7%). Meanwhile, compared with the
state-of-the-art tensor decomposed model-oriented hardware TIE, our proposed
FDHT-LSTM architecture achieves $2.5\times$, $1.46\times$ and $2.41\times$
increase in throughput, area efficiency and $2.41\times$ energy efficiency,
respectively on LSTM-Youtube workload. For LSTM-UCF workload, our proposed
design also outperforms TIE with $1.9\times$ higher throughput, $1.83\times$
higher energy efficiency and comparable area efficiency.

Figure 1: Architecture of tensor decomposition-based LSTM. (a) Prior TT-LSTM.
(b) The proposed FDHT-LSTM. Reproduced from [22].

The rest of this paper is organized as follows. Section 2 introduces the
background of tensor decomposition and the motivation of hierarchical Tucker
(HT)-based LSTM compression. The proposed fully decomposed HT-structured
(FDHT) LSTM algorithm is described in Section 3. Section 4 presents the
corresponding FDHT-LSTM hardware architecture. The evaluation performance at
both algorithm and hardware levels is reported in Section 5. Section 7 draws
the conclusions.

##  2 Background and Motivation

###  2.1 Preliminaries

####  2.1.1 Notation

In this paper, boldface lower-case, boldface capital and boldface calligraphic
letters represent vectors, matrices and higher-order tensors, respectively,
e.g. $\bm{x}$, $\bm{X}$ and $\bm{\mathcal{X}}$. Additionally, letters with
indices in bracket denote the entry, e.g., $\bm{x}(i)$, $\bm{X}(i,j)$,
$\bm{\mathcal{X}}(i_{1},\cdots,i_{d})$.

####  2.1.2 Tensor Contraction

Tensor contraction is the multiplication between two higher-order tensors
where more than one dimension matches. For example, given two tensors
$\bm{\mathcal{A}}\in\mathbb{R}^{M_{1}\times M_{2}\times C}$ and
$\bm{\mathcal{B}}\in\mathbb{R}^{C\times N_{1}\times N_{2}}$, where the 3rd
dimension of $\bm{\mathcal{A}}$ is identical to the 1st dimension of
$\bm{\mathcal{B}}$ with size $C$, the result of tensor contraction
$\bm{\mathcal{C}}=\bm{\mathcal{A}}\times_{1}^{3}\bm{\mathcal{B}}$, as a size-
$M_{1}\times M_{2}\times N_{1}\times N_{2}$ tensor, can be calculated as:

$\bm{\mathcal{C}}(i_{1},i_{2},j_{1},j_{2})=\sum_{k=1}^{C}\bm{\mathcal{A}}(i_{1},i_{2},k)\bm{\mathcal{B}}(k,j_{1},j_{2}).$
(1)

Figure 2: Standard HT decomposition example for a 4-order tensor. The
hierarchy is a binary tree with root $\\{1,2,3,4\\}$, where nodes are denoted
by rounded rectangles. $\\{1\\}$, $\\{2\\}$, $\\{3\\}$ $\\{4\\}$ are leaf
nodes, and $\\{1,2\\}$, $\\{3,4\\}$ are their parents. In HT format, only
colored leaf frames and transfer tensors are needed to store. Reproduced from
[22].

####  2.1.3 Hierarchical Tucker Decomposition

As a special case of tensor decomposition, hierarchical Tucker (HT) approach
recursively decomposes the original tensor into small tensor cores with
hierarchical levels from top to bottom in a binary tree. As illustrated in
Figure 2, the upper-level intermediate components can be factorized to lower-
level components. Here the intermediate components are referred as frames,
each of which corresponds to a unique node. Additionally, each node in the
binary tree is represented by a dimension set. In general, given a $d$-order
tensor $\bm{\mathcal{A}}\in\mathbb{R}^{N_{1}\times\cdots\times N_{d}}$, we can
build a binary tree to split the original dimension set $\\{1,\cdots,d\\}$ and
assign the index to each node. To be specific, the root node
$D=\\{1,2,\cdots,d\\}$ is associated with the root frame
$\bm{\mathcal{A}}=\bm{\mathcal{U}}_{D}$, $s\subsetneq D$ is the node
corresponding to frame $\bm{\mathcal{U}}_{s}$, and $s_{1},s_{2}\subsetneq s$
are the left and right child nodes of node $s$. Notice that if we define
$\mu_{s}=\min(s),\nu_{s}=\max(s)$, then each non-leaf frame
$\bm{\mathcal{U}}_{s}\in\mathbb{R}^{R_{s}\times N_{\mu_{s}}\times\cdots\times
N_{\nu_{s}}}$ can be recursively decomposed to frames
($\bm{\mathcal{U}}_{s_{1}}$ and $\bm{\mathcal{U}}_{s_{2}}$) and transfer
tensor $\bm{\mathcal{B}}_{s}\in\mathbb{R}^{R_{s}\times R_{s_{1}}\times
R_{s_{2}}}$ as

$\bm{\mathcal{U}}_{s}=\bm{\mathcal{B}}_{s}\times_{1}^{2}\bm{\mathcal{U}}_{s_{1}}\times_{1}^{2}\bm{\mathcal{U}}_{s_{2}},$
(2)

where $R_{s},R_{s_{1}},R_{s_{2}}$ are referred as hierarchical ranks. Overall,
by using hierarchical Tucker decomposition that recursively decomposes the
frames from top to bottom, we can use the combination of the small-size
2-order leaf frames and 3-order transfer tensors to store the original large-
size $d$-order tensor $\bm{\mathcal{X}}=\bm{\mathcal{U}}_{D}$.

###  2.2 Compressing LSTMs via Tensor Decomposition

In the recent years tensor decomposition has emerged as a very attractive
technique to compress very large LSTM models. In [20], tensor train (TT)
decomposition is applied to factorize the input-to-hidden layer of LSTMs and
GRUs for video recognition. The reported experimental results show that very
significant compression ratio can be achieved while maintaining high accuracy.
Motivated by this success, other advanced tensor decomposition approach, such
as block-term (BT) [21] decomposition and tensor ring (TR) decomposition [21],
are also used to build compact LSTM models (BT-LSTM and TR-LSTM) in video
processing tasks. In addition, tensor decomposed LSTMs also demonstrate high
performance in other sequential analysis and prediction tasks. For instance,
[23] proposes a TT-format LSTM to perform long-term forecasting in dynamic
systems.

###  2.3 Motivation of Compressing LSTM via HT Decomposition

In this paper we propose to adopt HT decomposition, a relatively little
noticed yet powerful tensor factorization approach for LSTM compression. Our
choice is motivated by two reasons. First, compared with its well-explored
counterparts (e.g. TT, TR and BT decomposition) in the model compression
field, HT decomposition, by it nature, can provide higher space complexity
reduction on the same-size tensor data with the same selected rank, thereby
implying that even smaller LSTM models can be constructed in the HT format.
Second, from the perspective of tensor theory, the inherent hierarchical
structure of HT also enables better weight sharing and hierarchical
representation from high-dimensional data, which are very important to
guarantee the representation capability of tensor decomposed LSTM models. As
we will report in Section 5, the HT structure-based LSTM models consistently
demonstrate superior performance than the existing compressed LSTMs using
other tensor decomposition methods with respect to both model accuracy and
compression ratio.

##  3 The Proposed FDHT-LSTM: Algorithm

###  3.1 Compact HT-structure Linear Layer

Consider the linear layer is the foundation of LSTM models, in this subsection
we first build a compact HT-structure layer that serves as the key component
of the proposed FDHT-LSTM model.

####  3.1.1 Tensorization

In the compact HT-layer, computations are performed with high-order tensor
contraction, thus we first need to transform all the original variables into
tensor format. In general, for an original linear layer with weight matrix
$\bm{W}\in\mathbb{R}^{O\times I}$ that linearly maps an input vector
$\bm{x}\in\mathbb{R}^{I}$ to an output vector $\bm{y}\in\mathbb{R}^{O}$ in HT-
based format, we first transform the weight matrix $\bm{W}$ to an $d$-order
tensor $\bm{\mathcal{W}}\in\mathbb{R}^{O_{1}\times\cdots\times O_{d}\times
I_{1}\times\cdots\times I_{d}}$. In addition, by defining
$I=\prod_{i=1}^{d}I_{i}$ and $O=\prod_{j=1}^{d}O_{j}$, we can also tensorize
the affiliated input vector $\bm{x}$ and output vector $\bm{y}$ to input and
output tensor $\bm{\mathcal{X}}\in\mathbb{R}^{I_{1}\times\cdots\times I_{d}}$
and $\bm{\mathcal{Y}}\in\mathbb{R}^{O_{1}\times\cdots\times O_{d}}$,
respectively.

####  3.1.2 Decomposing $\mathcal{W}$

After tensorizing the large-size weight matrix $\bm{W}$ as weight tensor
$\bm{\mathcal{W}}\in\mathbb{R}^{O_{1}\times\cdots\times O_{d}\times
I_{1}\times\cdots\times I_{d}}$, we can then utilize HT decomposition to
denote the original $\bm{W}$ with a set of small-size matrices and tensors.
With the definition of HT decomposition in Eq. (2), $\bm{\mathcal{W}}$ can be
factorized as

$\displaystyle\bm{\mathcal{W}}(o_{1},\cdots,o_{d},i_{1},\cdots,i_{d})=\sum_{k=1}^{R_{D}}\sum_{p=1}^{R_{D_{1}}}\sum_{q=1}^{R_{D_{2}}}\bm{\mathcal{B}}_{D}(k,p,q)\cdot$
(3)
$\displaystyle\bm{\mathcal{U}}_{D_{1}}(p,\varphi_{D_{1}}(\bm{o},\bm{i}))\bm{\mathcal{U}}_{D_{2}}(q,\varphi_{D_{2}}(\bm{o},\bm{i})),$

where $\varphi_{s}(\bm{o},\bm{i})$ is a mapping function which returns the
associate indices $\bm{o}=(o_{1},\cdots,o_{d})$ and
$\bm{i}=(i_{1},\cdots,i_{d})$ for a specified frame $\bm{\mathcal{U}}_{s}$
with a given node $s$ and $d$. For example, with $d=6$ and $s=\\{3,4\\}$, the
output of $\varphi_{s}(\bm{o},\bm{i})$ is $(o_{3},o_{4},i_{3},i_{4})$.
Additionally, $\bm{\mathcal{U}}_{D_{1}}$ and $\bm{\mathcal{U}}_{D_{2}}$ can be
further recursively decomposed as

$\displaystyle(\bm{\mathcal{U}}_{s})(k,\varphi_{s}(\bm{o},\bm{i}))=\sum_{p=1}^{R_{s_{1}}}\sum_{q=1}^{R_{s_{2}}}(\bm{\mathcal{B}}_{s})(k,p,q)$
(4)
$\displaystyle\cdot(\bm{\mathcal{U}}_{s_{1}})(p,\varphi_{s_{1}}(\bm{o},\bm{i}))(\bm{\mathcal{U}}_{s_{2}})(q,\varphi_{s_{1}}(\bm{o},\bm{i})),$

where $D=\\{1,2,\cdots,d\\}$, $D_{1}=\\{1,\cdots,\lfloor d/2\rfloor\\}$ and
$D_{2}=\\{\lceil d/2\rceil,\cdots,d\\}$ are the left and right child nodes of
the root node $D$.

####  3.1.3 HT-Structure Layer

As the original weight matrix is decomposed to HT-format, the HT-structured
linear layer can be readily constructed via HT-format matrix-vector
multiplication instead of recovering back to original format. To be specific,
the forward pass of an HT-layer is computed as:

$\displaystyle\bm{\mathcal{Y}}(\bm{o})=\sum_{\bm{i}}\sum_{k=1}^{R_{D}}\sum_{p=1}^{R_{D_{1}}}\sum_{q=1}^{R_{D_{2}}}(\bm{\mathcal{B}}_{D})(k,p,q)\cdot$
(5)
$\displaystyle\bm{\mathcal{U}}_{D_{1}}(p,\varphi_{D_{1}}(\bm{o},\bm{i}))\bm{\mathcal{U}}_{D_{2}}(q,\varphi_{D_{2}}(\bm{o},\bm{i}))\bm{\mathcal{X}}(\bm{i}).$

After the above computation, we need to reshape the output tensor
$\bm{\mathcal{Y}}$ back to the original vector-format $\bm{y}$. Hence, we can
simplify the entire computing process as a single HT-structure layer:

$\bm{y}=HTL(\bm{W},\bm{x}).$ (6)

As depicted in Figure 3, the consecutive arrows denote the computation flow in
the computation process of an HT-structure layer. Specifically, the input
vector $\bm{x}$ is first tensorized, and then it is contracted with the HT-
format weight tensor in a hierarchical way. Finally, the obtained output
tensor is reshaped back to output vector $\bm{y}$.

Figure 3: Computation flow in the HT-layer for $d$=4. Arrows represent the
directions of tensor contractions among the decomposed tensors. Here leaf
frames $\\{\bm{\mathcal{U}}\\}_{i=1}^{4}$ are 3-order tensors since the weight
tensor is tensorized from a matrix. Reproduced from [22].

####  3.1.4 Benefits on Low Cost

With the proposed HT-layer, one most important benefit is huge complexity
reduction. Table I shows the theoretical space complexity in comparison with
other tensor factorized linear layer as well as the uncompressed one.
Considering tensor rank $R$ is generally smaller than $I$ or $O$, HT-structure
layer can provide the lowest space complexity as shown in this table.
Additionally, to verity the practical ability of parameters reduction, we plot
curves of number of parameters with respect to $R$ in different tensor
decomposition formats. As shown in Figure 4, to store size-$57,600\times 256$
weight matrix used in [20, 21, 17], our HT-structure layer requires the fewest
number of parameters with the same rank settings among all different types of
tensor factorization methods. As will be shown in Section 5, such advantages
will further be verified via empirical experimental results across various
popular video recognition datasets.

TABLE I: Space complexities for different tensor-format linear layers. Here $R=\max_{s\subsetneq D}R_{s}$, $I^{\prime}=\max_{k\in D}I_{k}$, $O^{\prime}=\max_{k\in D}O_{k}$. Reproduced from [22]. Compressed Linear Layer | Space Complexity
---|---
Uncompressed | $\mathcal{O}(I^{\prime}O^{\prime})$
Tensor Train (TT)-structure | $\mathcal{O}(dI^{\prime}O^{\prime}R^{2})$
Tensor Ring (TR)-structure | $\mathcal{O}(dI^{\prime}O^{\prime}R^{2})$
Block-Term (BT)-structure | $\mathcal{O}(dI^{\prime}O^{\prime}R+R^{d})$
Hierarchical Tucker (HT)-structure | $\mathcal{O}(dI^{\prime}O^{\prime}R+dR^{3})$

###  3.2 Fully Decomposing the HT-structure LSTM

####  3.2.1 Challenges of Fully Compressing LSTM

Recall that a typical LSTM model consists of input-to-hidden layer and hidden-
to-hidden layer. Built on the top of the underlying HT structure for each
component layer, performing fully decomposition on the entire LSTM models can
evidently bring further performance improvement. To that end, one naive way is
to simply compress each layer with HT decomposition. In other words, all
weight matrices in the LSTM model are independently decomposed to HT format.
Although this strategy can indeed bring full compression on the entire LSTM
network, such straightforward layer-wise compression strategy suffers from
huge accuracy drop. Figure 5 shows the experimental results using the layer-
wise compression with HT decomposition. Here more than $2\%$ accuracy drop can
be observed as compared to the conventional input-to-hidden-only compression
adopted in [20, 21, 17]. Overall, fully compressing the entire LSTM model
using HT decomposition without accuracy drop is challenging and non-trivial.

Figure 4: Curves for number of parameters in different tensor-format layer
with respect to tensor rank $R$. Here we use settings in the related works,
i.e., $d=5$, $(I_{1},\cdots,I_{5})=(8,10,10,9,8)$,
$(O_{1},\cdots,O_{5})=(4,4,2,4,2)$. Reproduced from [22]. Figure 5: Comparison
between naive layer-wise compression and our proposed full decomposition in
terms of number of parameters and test accuracy on UCF11 dataset. Reproduced
from [22].

####  3.2.2 Proposed FDHT-LSTM

To fully compress LSTM models without performance degradation, we propose a
homogeneous compression method, namely fully decomposed HT (FDHT)
decomposition, to obtain highly-compact LSTM models. Figure 1 illustrates the
key idea of the proposed full decomposition. To be specific, in order to
maximally leverage the linear correlation across all weight matrices in the
entire LSTM model, we first concatenate all the weight matrices as a single
huge matrix such that the entire model can be considered as a “mega” linear
layer. Then, at each time step, all the intermediate results only need one-
time multiplication in the forward propagation process:

$\displaystyle\begin{bmatrix}\hat{\bm{f}}[t]\\\ \hat{\bm{u}}[t]\\\
\hat{\bm{c}}[t]\\\ \hat{\bm{o}}[t]\end{bmatrix}$
$\displaystyle=\begin{bmatrix}\bm{W}_{f}&\bm{V}_{f}\\\
\bm{W}_{u}&\bm{V}_{u}\\\ \bm{W}_{c}&\bm{V}_{c}\\\
\bm{W}_{o}&\bm{V}_{o}\end{bmatrix}\begin{bmatrix}\bm{x}[t]\\\
\bm{h}[t-1]\end{bmatrix}$ (7) $\displaystyle=\bm{W}\bm{I}[t].$

With the above interpretation, we can directly impose the desired HT structure
on the integrated single huge matrix in an LSTM model. To be specific,
following the scheme described in Section 3.1, we can tensorize and decompose
the entire LSTM models to the HT format, and then perform forward pass as:

$\begin{bmatrix}\hat{\bm{f}}[t]\\\ \hat{\bm{u}}[t]\\\ \hat{\bm{c}}[t]\\\
\hat{\bm{o}}[t]\end{bmatrix}=\bm{Z}[t]=HTL(\bm{W},\bm{I}[t]).$ (8)

Correspondingly, the outputs of the FDHT-LSTM can be calculated as follows:

$\displaystyle\bm{f}[t]=$ $\displaystyle\sigma(\hat{\bm{f}}[t]+\bm{b}_{f})$
(9) $\displaystyle\bm{u}[t]=$
$\displaystyle\sigma(\hat{\bm{u}}[t]+\bm{b}_{u})$ $\displaystyle\bm{c}[t]=$
$\displaystyle\bm{f}[t]\odot\bm{c}[t-1]+\bm{u}[t]\odot\tanh(\hat{\bm{c}}[t]+\bm{b}_{c})$
$\displaystyle\bm{o}[t]=$ $\displaystyle\sigma(\hat{\bm{o}}[t]+\bm{b}_{o})$
$\displaystyle\bm{h}[t]=$ $\displaystyle\bm{o}[t]\odot\tanh(\bm{c}[t]),$

where $\sigma$, $\tanh$ and $\odot$ are the sigmoid function, hyperbolic
function and element-wise product, respectively.

####  3.2.3 HT-based Gradient Calculation

Notice that in order to ensure that the desired FDHT-LSTM model can be
properly obtained, the backward propagation in the training process should
also be reformulated to HT-based format. In general, given an HT-structure
linear layer $\bm{\mathcal{W}}=\bm{\mathcal{U}}_{D}$ and
$\frac{\partial\bm{\mathcal{Y}}}{\partial\bm{\mathcal{U}}_{D}}=\bm{\mathcal{X}}$,
and defining $s$, $F(s)$ and $B(s)$ as non-root node, parent and sibling nodes
of node $s$ in the binary tree, respectively, the partial derivative of the
output tensor with respect to frames can be recursively calculated as:

$\displaystyle\frac{\partial\bm{\mathcal{Y}}}{\partial\bm{\mathcal{U}}_{s}}=$
$\displaystyle\bm{\mathcal{B}}_{F(s)}\times_{1}^{3}\bm{\mathcal{U}}_{B(s)}$
(10)
$\displaystyle\times_{\begin{subarray}{c}1,\cdots,\nu_{B(s)}-\mu_{B(s)}+2,\nu_{F(s)}-\mu_{F(s)}+3,\\\
\cdots,\nu_{F(s)}-\mu_{F(s)}+\nu_{B(s)}-\mu_{B(s)}+3\end{subarray}}^{1,3,\cdots,2\nu_{B(s)}-2\mu_{B(s)}+4}\frac{\partial\bm{\mathcal{Y}}}{\partial\bm{\mathcal{U}}_{F(s)}},$

where $\mu_{s}=\min(s),\nu_{s}=\max(s)$. Notice that when $F(s)$ is equal to
node $D$, the above recursive procedure terminates.

Furthermore, with Eq. (10) the gradients for leaf frames and transfer tensors
can also be recursively calculated as:

$\displaystyle\frac{\partial L}{\partial\bm{\mathcal{U}}_{s}}=$
$\displaystyle\frac{\partial\bm{\mathcal{Y}}}{\partial\bm{\mathcal{U}}_{s}}\times_{1,\cdots,\mu_{s}-1,\nu_{s}+1,\cdots,d}^{\nu_{s}-\mu_{s}+3,\cdots,d+1}\frac{\partial
L}{\partial\bm{\mathcal{Y}}},$ (11) $\displaystyle\frac{\partial
L}{\partial\bm{\mathcal{B}}_{s}}=$
$\displaystyle\frac{\bm{\mathcal{Y}}}{\partial\bm{\mathcal{U}}_{s}}\times_{2,\cdots,\nu_{s_{1}}-\mu_{s_{1}}+2}^{2,\cdots,\nu_{s_{1}}-\mu_{s_{1}}+2}\bm{\mathcal{U}}_{s_{1}}$
$\displaystyle\times_{2,\cdots,\nu_{s_{2}}-\mu_{s_{2}}+2}^{3,\cdots,\nu_{s_{2}}-\mu_{s_{2}}+3}\bm{\mathcal{U}}_{s_{2}}\times_{1,\cdots,d}^{4,\cdots,d+3}\frac{\partial
L}{\partial\bm{\mathcal{Y}}}.$ (12)

In general, our proposed “integrate-then-decompose” compression strategy can
bring significant performance benefit for the FDHT-LSTM models. As shown in
Figure 5, compared with the naive layer-wise compression strategy that suffers
significant accuracy degradation, our proposed approach ensures high accuracy
without performance loss. Furthermore, compared to the input-to-hidden-only
compression counterpart, our FDHT-LSTM enjoys much smaller model size. More
detailed and comprehensive empirical evaluation results across different
datasets will be reported in Section 5.

##  4 The Proposed FDHT-LSTM: Hardware Architecture

Based on the above described FDHT-LSTM structure and algorithm, in this
section we further develop the corresponding hardware architecture that can
fully leverage the algorithmic benefits to support the efficient execution of
the FDHT-LSTM model.

###  4.1 Analysis of Computation Flow: A 2-D Matrix Perspective

Figure 6: Computation flow of high-order FDHT-LSTM from the perspective of 2-D
matrix multiplication and transformation.

As described in Section 3, the overall computation scheme of FDHT-LSTM can be
interpreted as a series of tensor contraction along a binary tree (see Figure
3). Consider tensor contraction (Eq. 1) is essentially a high-dimensional
operation; while the underlying arithmetic and memory hardware components can
only support 2-D operation, a re-analysis of the computational flow, from 2-D
matrix perspective, is very necessary and desired.

Figure 6 illustrates the overall computational flow of the proposed FDHT-LSTM
in the 2-D matrix format. Here $\bm{B}^{{}^{\prime}}_{1234}$,
$\bm{B}^{{}^{\prime}}_{12}$, $\bm{B}^{{}^{\prime}}_{34}$,
$\bm{U}^{{}^{\prime}}_{1}$, $\bm{U}^{{}^{\prime}}_{2}$,
$\bm{U}^{{}^{\prime}}_{3}$ and $\bm{U}^{{}^{\prime}}_{4}$ are the weight
matrices of the flatten high-order tensors, $\bm{X^{{}^{\prime}}}$ is the
reshaped input matrix data, and $\bm{T}$ is the intermediate result. Notice
that here $\bm{T}$ cannot be multiplied with the next weight matrix directly
due to the mismatch of matrix dimension. Therefore, $\bm{T}$ needs to be first
transformed to $\bm{T}^{{}^{\prime}}$ to ensure functional validity.

As shown in Figure 6, a total of six matrix transformations (from $\bm{T}$ to
$\bm{T^{{}^{\prime}}}$) is required in the entire computation flow of FDHT-
LSTM. After further analysis, we can obtain two important observations: First,
these transformations are not simple matrix reshape but essentially
permutation that requires complex matrix operations such as row/column merge
and concatenation etc.; Second the six transformations can be categorized to
three types. To be specific, Type-I (Transformation 1 and 4) is to convert
$\bm{T}\in\mathbb{R}^{A\times(B_{1}\times B_{2})}$ to
$\bm{T^{{}^{\prime}}}\in\mathbb{R}^{(A\times B_{1})\times B_{2}}$, and the
mapping principle for this transformation is as:

$\displaystyle\begin{split}\bm{T}&(m,n)\Rightarrow\bm{T^{{}^{\prime}}}(p,q),\end{split}$
(13)

where $m=a$, $n=(b_{1}-1)\times B_{2}+b_{2}$, $p=(a-1)\times B_{1}+b_{1}$,
$q=b_{2}$, $a=1,2,...,A$, $b_{1}=1,2,...,B_{1}$ and $b_{2}=1,2,...,B_{2}$.

Type-II (Transformation 2 and 5) is to convert
$\bm{T}\in\mathbb{R}^{(A_{1}\times A_{2}\times A_{3})\times(B_{1}\times
B_{2})}$ to $\bm{T^{{}^{\prime}}}\in\mathbb{R}^{(A_{1}\times A_{3}\times
B_{1})\times(A_{2}\times B_{2})}$, and the mapping principle for this
transformation is as:

$\displaystyle\begin{split}\bm{T}&(m,n)\Rightarrow\bm{T^{{}^{\prime}}}(p,q),\end{split}$
(14)

where $m=(a_{1}-1)\times A_{2}\times A_{3}+(a_{2}-1)\times A_{3}+a_{3}$,
$n=(b_{1}-1)\times B_{2}+b_{2}$, $p=(a_{1}-1)\times A_{3}\times
B_{1}+(a_{3}-1)\times B_{1}+b_{1}$, $q=(a_{2}-1)\times b_{2}$,
$a_{1}=1,2,...,A_{1}$, $a_{2}=1,2,...,A_{2}$, $a_{3}=1,2,...,A_{3}$,
$b_{1}=1,2,...,B_{1}$ and $b_{2}=1,2,...,B_{2}$.

Type-III (Transformation 3 and 6) is to convert
$\bm{T}\in\mathbb{R}^{(A_{1}\times A_{2}\times A_{3})\times B}$ to
$\bm{T^{{}^{\prime}}}\in\mathbb{R}^{(A_{1}\times A_{3}\times B)\times A_{2}}$,
and the mapping principle for this transformation is as:

$\displaystyle\begin{split}\bm{T}&(m,n)\Rightarrow\bm{T^{{}^{\prime}}}(p,q),\end{split}$
(15)

where $m=(a_{1}-1)\times A_{2}\times A_{3}+(a_{2}-1)\times A_{3}+a_{3}$,
$n=b$, $p=(a-1)\times A_{3}\times B+(a_{3}-1)\times B+b$, $q=a_{2}$,
$a_{1}=1,2,...,A_{1}$, $a_{2}=1,2,...,A_{2}$, $a_{3}=1,2,...,A_{3}$ and
$b=1,2,...,B$.

Figure 7: Data access in conventional SRAM architecture.

###  4.2 Access Conflict on Conventional SRAM Architecture

Without losing the generality, we use the basic matrix transformation
$\bm{M}\in\mathbb{R}^{(A_{1}\times A_{2})\times(B_{1}\times
B_{2})}\Rightarrow\bm{M}^{{}^{\prime}}\in\mathbb{R}^{(A_{1}\times
B_{1})\times(A_{2}\times B_{2})}$ to elaborate the challenges incurred by the
complicated transformation above, where $A_{1}=A_{2}=B_{1}=B_{2}=2$. As shown
in Figure 7, when using conventional SRAM architecture, in the writing phase,
each row of the matrix is written to one address (one row) in this SRAM in one
cycle. For example, data $01,02,03$ and $04$ is written to the address $01$ of
the SRAM. However, the access conflict happens in the read phase. When reading
the SRAM, the desired data is $01,02$ and $04,05$. These four data are located
at two different addresses $01$ and $02$, and these two addresses can not be
accessed simultaneously, thereby causing access conflict.

The basic solution to the access conflict problem is relocating one address
multiple times (as shown in Figure 7). More specifically, in cycle 1, address
$01$ is located and data $01,02,03,04$ is read out from the SRAM. Data $01,02$
is saved in the register file, while $03,04$ is the undesired data and
discarded. Similarly, in cycle 2, address $02$ is located; $05,06,07,08$ is
read out from SRAM, and $05,06$ is stored in the register file and $07,08$ is
discarded. In the assemble unit, $01,02$ and $05,06$ form the data array
$01,02,05,06$. However, such a solution suffers two main limitations:

1\. Additional memory overhead is needed to support matrix transformation.
Since the matrix transformation cannot be executed in an on-the-fly way when
using conventional memory architecture, an additional register file is
necessary to save the desired data. Thus, this solution significantly
increases the cost of the register file, reducing the area and energy
performance of the overall hardware.

2\. The write and reading scheme becomes much more complicated because of
frequent re-locating of the memory address. Also, notice that the example we
raise here is only a simple case, and in practical workload, the matrix
transformation is very large and complicated. Hence it is not feasible to read
all the data to register files and assemble them when using conventional SRAM
architecture.

###  4.3 2-D SRAM Array

In order to implement the desired matrix transformation, which is essentially
the 2-D mapping of high-order tensor transpose, we propose to read and write
the data in the 2-D SRAM array. In other words, we aim to increase the
dimension of the SRAM from physical 2-D to ”virtually” 4-D, to better support
the operation in the high-order tensor operation. As illustrated in Figure 8,
our proposed SRAM array physically consists of $G$ banks of SRAM, where each
of them has width of $W$ and depth of $M$. Based on such fixed arrangement,
the controller will logically partition the SRAM array to more fine-grained
format when supporting different workloads with various matrix sizes. To be
specific, according to the demand of the workload, each SRAM bank is logically
partitioned into $N$ segments along the depth dimension, and each segment has
a depth of $D$ where $M=N\times D$. Notice that here $G$, $M$ and $W$ are the
fixed dimensions of the entire SRAM module, and $D$ and $N$ are reconfigurable
and determined by the main controller, thereby providing the controller
sufficient capability to adjust the granularity of the memory operation for
different workloads.

To achieve this goal without introducing any memory access conflict, delicate
design of memory read/write scheme is required and will be described next.

Figure 8: The physical SRAM bank and logically partitioned 2-D SRAM array.
Here $W$ and $M$ are the width and depth of one SRAM bank, respectively. $D$
is the depth of one segment of SRAM bank, where $M=D\times N$. Here $D$ can be
reconfigured according to workload.

###  4.4 Matrix Transformation on 2-D SRAM Array

Based on the 2-D SRAM array, the corresponding conflict-free memory access
scheme can be developed. Without loss of generality, in this subsection we
first study the read and write schemes for a basic transformation
$\bm{M}\in\mathbb{R}^{(A_{1}\times A_{2})\times(B_{1}\times
B_{2})}\Rightarrow\bm{M}^{{}^{\prime}}\in\mathbb{R}^{(A_{1}\times
B_{1})\times(A_{2}\times B_{2})}$. The obtained outcome will further serve as
the foundation of implementing the specific matrix transformations (Eq.13 –
Eq.15) used in FDHT-LSTM.

To facilitate the notation, as illustrated in Figure 8, the bank and segment
indices of the 2-D SRAM array are denoted as x and y, respectively, where
$x=1,2,...,G$ and $y=1,2,...,N$. In addition, z and k are the row and depth
indices in each SRAM segment, respectively, where $z=1,2,...,W$ and
$k=1,2,...,D$. Following such notation, one data entry in one SRAM segment,
one row in one SRAM segment, one SRAM segment, and one SRAM bank can be
represented as $S(x,y,k,z)$, $S(x,y,k,:)$, $S(x,y,:,:)$ and $S(x,:,:,:)$,
respectively. Also, we use $\bm{M}(a,:)$ and $\bm{M}(:,b)$ to represent the
$a$-th row vector and the $b$-th column vector of $\bm{M}$, respectively. In
addition, $\bm{M}(a,b:c)$ and $\bm{M}(a:b,c)$ denote the part of $\bm{M}(c,:)$
and $\bm{M}(:,c)$ that are with column and row indices from $a$ to $b$.

Figure 9: The proposed write and read scheme to 2-D SRAM array for a basic
matrix transformation.

####  4.4.1 Memory Write Scheme

Next we describe memory write scheme for the example transformation
$\bm{M}\in\mathbb{R}^{(A_{1}\times A_{2})\times(B_{1}\times
B_{2})}\Rightarrow\bm{M}^{{}^{\prime}}\in\mathbb{R}^{(A_{1}\times
B_{1})\times(A_{2}\times B_{2})}$. Here $D$, as the granularity of 2-D SRAM
array, is configured as $B_{2}$ to support the operation for this example.
Initially, $\bm{M}(1,1:B_{2})$ is written along the $z$ dimension into
$S(1,1,1,:)$ as the first row of SRAM segment $S(1,1)$. Then, if we interpret
the $B_{2}$ data entries in one row as a group, and one row of $\bm{M}$
contains $B_{1}$ groups, these $B_{1}$ data groups are scheduled to be
sequentially written into $S(1,1,:,:)$ in a row-wise way. For instance, the
data entries of $\bm{M}(1,B_{2}+1:2\times B_{2})$ are written into
$S(1,1,2,:)$, and the data entries of $\bm{M}(1,(B_{1}-1)\times
B_{2}+1:B_{1}\times B_{2})$ are written into $S(1,1,B_{1},:)$. After that, a
row of data entries of $\bm{M}$ is stored in one segment of 2-D SRAM array.
Then, $A_{2}$ rows of $\bm{M}$ are written into $A_{2}$ SRAM segments along
the x dimension, e.g., $S(1,2,:,:)$ stores the data of $\bm{M}(2,:)$ and
$S(1,S_{2},:,:)$ stores the data of $\bm{M}(A_{2},:)$. By using this way,
$\bm{M}(1:A_{2},:)$ is written into one SRAM bank $S(1,:,:,:)$. Writing the
rest of the data in $\bm{M}$ to SRAM follows the similar way, e.g.,
$S(2,:,:,:)$ stores the data of $\bm{M}(A_{2}+1:2\times A_{2},:)$, and
$S(A_{1},:,:,:)$ stores the data of $\bm{M}((A_{1}-1)\times
A_{2}+1:A_{1}\times A_{2},:)$. Figure 9 shows the detailed scheme for an
example matrix transformation, where $G=2$, $M=4$, $W=2$, and
$A_{1}=A_{2}=B_{1}=B_{2}=2$. It is seen that with such writing scheme the 2-D
SRAM array stores the re-arranged data, which will be ready to realize the
desired transformation after proper reading scheme.

####  4.4.2 Memory Read Scheme

In the reading phase, the data in $G$ banks are read simultaneously in a row-
wise way. To be specific, the data in $S(:,1,:,k)$ are first read in the first
cycle, and then they are sent to assemble unit to form
$\bm{M^{{}^{\prime}}}(1,:)$. In the following cycles, such operation is
repeated row by row along the k and y dimensions. Figure 9 shows the details
of reading scheme for the example matrix transformation.

In general, the proposed 2-D SRAM with a corresponding write and read scheme
can solve the challenges mentioned above:

1\. The matrix transformation can be executed in an on-the-fly way. In the
reading phase, all the to-be-read data is the desired data and only a small
register file is needed to save the data. In our FDHT-LSTM, the size of the
register file is only $448Bytes$.

2\. Since in the writing and reading phase, the desired memory address
increases sequentially, it is very easy to design and implement the address
controller. Also, because the 2-D SRAM architecture is logically split, it is
very easy to be reconfigured and very flexible for different matrix
transformations with various shapes.

###  4.5 Realizing Matrix Transformation for FDHT-LSTM

Recall that the read/write access schemes described in Section 4.4 are
designed for a basic matrix transformation that is different from the ones
described in Eq. 13 – Eq. 15. In this subsection we leverage the approach
developed in Section 4.4 to further realize the three matrix transformations
specifically adopted in the computation flow of FDHT-LSTM.

####  4.5.1 Memory Access Scheme for Type-I Transformation

Figure 10: Example memory write and read schemes for Transformation 1.

In this case $D$ is set as 1. During the memory writing phase, consider
$B_{2}$ data as a group, then the data in each group are written to one row of
SRAM segment along the z dimension. Since each row of $\bm{T}$ contains
$B_{1}$ data groups, they are written into one row of SRAM array along the x
dimension sequentially. In other words, one row of SRAM array $S(:,1,:,:)$
stores one row of $\bm{T}$ as $\bm{T}(1,:)$. Following this strategy, the $A$
rows of $\bm{T}$ are written into the SRAM array row by row along the y
dimension in a sequential way. Then during the memory reading phase, consider
each row of SRAM array $S(:,1,:,:)$ has $B_{1}\times B_{2}$ data; hence the
assemble unit will partition one row data to form a data array with $B_{1}$
columns and $B_{2}$ rows. By repeating such reading operation along the $y$
dimension, $\bm{T^{{}^{\prime}}}$ can be read out and sent to PE array for
further processing. Figure 10 illustrates the details of memory write and read
schemes to an example Transformation 1, where $A=I_{1}\times I_{2}\times
I_{3}=3$, $B_{1}=O_{4}=3$ and $B_{2}=R_{4}=2$.

####  4.5.2 Memory Access Scheme for Type-II Transformation

Figure 11: Example memory write and read schemes for Transformation 2.

In this case, the transformation can be rewritten as
$\bm{T}\in\mathbb{R}^{(A_{1}\times A_{2})\times((A_{3}\times B_{1})\times
B_{2})}\Rightarrow\bm{T^{{}^{\prime}}}\in\mathbb{R}^{(A_{1}\times(A_{3}\times
B_{1}))\times(A_{2}\times B_{2})}$, and then it becomes the basic format
example described in Section 4.4 with $D=A_{3}\times B_{1}$. Then, during the
writing phase, we rearrange the data of $\bm{T}$ along $B_{2}$, $A_{3}\times
B_{1}$, $A_{2}$ and $A_{1}$ dimensions on the z, k, x and y dimensions of 2-D
SRAM array, respectively. In the memory reading phase, $A_{2}\times B_{2}$
data entries in the SRAM array are read simultaneously and combined as one row
of $\bm{T^{{}^{\prime}}}$. In general, by reading along the $k$ and $y$
dimensions of 2-D SRAM array, the row data of $\bm{T^{{}^{\prime}}}$ as
$A_{1}\times A_{3}\times B_{1}$ entries can be obtained. Figure 11 illustrates
the details of memory access scheme to an example Transformation 2, where
$B_{2}=R_{3}=2$, $A_{3}\times B_{1}=O_{4}\times R_{34}=3$, $A_{2}=I_{3}=3$ and
$A_{1}=I_{1}\times I_{2}=2$.

####  4.5.3 Memory Access Scheme for Type-III Transformation

Figure 12: Example memory write and read schemes for Transformation 3.

In the writing phase of this case, consider each row of $\bm{T}$ as a group,
and the data in each group is written along the dimension z in each SRAM
segment, e.g., $S(1,1,1,:)$ stores the data in $\bm{T}(1,:)$. Consequently, in
order to keep $A_{3}$ in the column index of $\bm{T^{{}^{\prime}}}$, $D$
should be set as $A_{3}$ to ensure that each $A_{3}$ rows of $\bm{T}$ are
written along the $k$ dimension in each SRAM segment, e.g.,
$\bm{T}(1:A_{3},:)$ is stored in $S(1,1,:,:)$. In addition, the data of
$\bm{T}$ along $A_{2}$ and $A_{1}$ dimensions are re-arranged on the x and y
dimensions of 2-D SRAM array, respectively, to guarantee $A_{2}$ appear in the
row index of $\bm{T}^{{}^{\prime}}$. During the reading phase, $A_{2}\times B$
data entries in the SRAM array are read simultaneously, and they are reshaped
into an array with $A_{2}$ columns and $B$ rows in the assemble units. Figure
12 illustrates the details of memory read and write schemes of the
Transformation 3, where $B=O_{3}=2$, $A_{3}=O_{4}\times R_{34}=2$, and
$A_{1}=I_{1}=1$.

###  4.6 Overall Architecture

Figure 13: The overall hardware architecture.

Based on the proposed 2-D SRAM array and the corresponding read/write access
schemes, the hardware architecture that supports the execution of FDHT-LSTM
models can be developed. As illustrated in Figure 13, the datapath of the
accelerator is a PE Array that consists of multiply-accumulators (MACs) to
perform matrix multiplication. Notice that because FDHT-LSTM operates on the
the dense matrices rather than the sparse ones, the multiplier utilization is
very high. In addition, the weight of the factorized tensor cores and the
intermediate results are stored in the weight SRAM and working SRAM,
respectively. Notice that here only working SRAM is in the format of 2-D array
since the desired matrix transformation is only involved with the intermediate
results. The read/write access schemes described in Section 4.5 are mapped
into an individual address generator, which receives the configuration
information from main controller for different workloads.

##  5 Evaluation

###  5.1 Performance of FDHT-LSTM Model

We evaluate the performance of our proposed FDHT-LSTM approach on various
video recognition tasks and compare it with other state-of-the-art tensor
decomposition-based LSTM compression methods. In particular, in order to
verify the benefits of our full decomposition strategy and fairly compare the
HT decomposition with other tensor decomposition for LSTM compression, we also
conduct experiments on only performing HT decomposition on the input-to-hidden
layer of LSTM models.

####  5.1.1 Experimental Settings

For fair comparison, we set $d$, the order of the target tensors of the
tensorized input, output and weights, identical to the settings in prior works
[20] and [21]. For the tensor ranks, we also follow the prior works to select
$r$ towards reaching a good balance between compression ratio and model
accuracy. For the setting of other hyper-parameters, dropout rate and batch
size are set as 0.25 and 16, respectively, and ADAM is selected as the
underlying training optimizer with weight decay of 0.001.

TABLE II: Compression settings for HT-layer. Here because YTC dataset is smaller than UCF11 dataset, we choose smaller $R$ (non-leaf) for LSTM-YTC to achieve higher compression ratio with still preserving high accuracy. Model | LSTM-UCF11 | LSTM-YTC
---|---|---
Size | $(57600,256)$ | $(57600,256)$
$d$ | $4$ | $4$
$I$ | $[16,16,16,15]$ | $[16,16,16,15]$
$O$ | $[4,4,4,4]$ | $[4,4,4,4]$
$R$ (Leaf) | $14$ | $14$
$R$ (Non-leaf) | $12$ | $11$
Compr. Ratio | 6,726$\times$ | 7,117$\times$

####  5.1.2 Results on UCF11 Dataset

The UCF11 dataset [16] is a set of 1,600 human action video clips which can be
categorized as 11 classes. The resolution of each clip is $320\times 240$. We
follow the prior works [21] [17] to pre-process the video data. To be
specific, the resolution of each video clip is downgraded to $160\times 120$,
then we randomly sample 6 frames from each clip and combine them as a single
sequential data point. In the original uncompressed LSTM model, there are 4
input-to-hidden layers and 4 hidden-to-hidden layers. We also follow the prior
works to set the input size and the number of hidden states in the baseline
model as $160\times 120\times 3=57,600$ and 256, respectively. For the
compressed FDHT-LSTM model, the concatenated vector $\bm{I}$ can be reshaped
as a tensor of shape as $16\times 16\times 16\times 15$, and the hidden state
is reshaped to a tensor of shape as $4\times 4\times 4\times 4$. In addition,
all leaf and non-leaf ranks are set as 14 and 12, respectively.

TABLE III: Compression performance of different tensor decomposed-based LSTM models on UCF11 dataset. Model | Number of Parameters | Top-1
---|---|---
Input-to-Hidden | Overall |  Acc. (%)
LSTM | 58.98M | 59.24M | 69.7
TT-LSTM [20] | 3,360 | 265.50K | 79.6
(ICML’17) | (223$\times$)
BT-LSTM [21] | 3,387 | 265.53K | 85.3
(CVPR’18) | (223$\times$)
TR-LSTM [17] | 1,725 | 263.87K | 86.9
(AAAI’19) | (225$\times$)
HT-LSTM | 1,245 | 263.39K | 87.2
(Ours) | (225$\times$)
FDHT-LSTM | N/A | 8,808 | 87.5
(Ours) | (6,726$\times$)

The compression performance of the proposed FDHT-LSTM and other state-of-the-
art tensor decomposed LSTM are summarized in Table III. It is seen that
compared with the original uncompressed LSTM model with 59 million parameters,
our proposed FDHT-LSTM only needs 8,808 parameters while resulting in a 17.8%
accuracy increase. Compared with other tensor decomposition-based LSTM
compression methods, i.e., TT-LSTM, BT-LSTM and TR-LSTM, FDHT-LSTM enjoys
order-of-magnitude reduction in model parameters and at least 0.6% higher
accuracy.

####  5.1.3 Results on Youtube Celebrities Face Dataset

Youtube celebrities face dataset [13] contains 1,910 video clips with
different resolutions which are categorized as 47 labels. Similar to the case
of UCF11 dataset, we also follow the prior works [21] [17] for data pre-
processing. To be specific, we scale down all the video clips to $160\times
120$ and randomly sample 6 frames to form a single sequence. Other
experimental settings are also consistent with the case of UCF11 dataset
except that all the non-leaf ranks are set as 11.

Table IV summarizes the compression performance on this dataset for different
tensor decomposed LSTMs. It is seen that compared to the baseline LSTM model
our FDHT-LSTM can provide as high as 7,117$\times$ compression ratio with test
accuracy increase. Compared with the state-of-the-art TT-LSTM, our FDHT-LSTM
still enjoys 6,894$\times$ smaller model size and 12.7% higher accuracy.

TABLE IV: Compression performance of different tensor decomposed-based LSTM models on Youtube celebrities face dataset. Model | Number of Parameters | Top-1
---|---|---
Input-to-Hidden | Overall |  Acc. (%)
LSTM | 58.98M | 59.24M | 33.2
TT-LSTM [20] | 3,392 | 265.54K | 75.5
(ICML’17) | (223$\times$)
HT-LSTM | 810 | 262.95K | 88.1
(Ours) | (225$\times$)
FDHT-LSTM | N/A | 8,324 | 88.2
(Ours) | (7,117$\times$)

Besides, we also compare the performance of FDHT-LSTM with several other
tensor decomposition-free models for video recognition, and the results are
summarized in Table V. It is seen that our FDHT-LSTM can outperform the state-
of-the-art work [14] with 3.6% higher accuracy and much smaller model size.

TABLE V: Experimental comparison among FDHT-LSTM and other non-tensor decomposition models, e.g., DML-PV [5], VGGFACE + RRNN [14] and VGG16-GCR [15], on Youtube celebrities face dataset. Model |  Number of Parameters | Top-1
---|---|---
Acc. (%)
DML-PV | 220K | 82.8
VGGFACE + RRNN |  $\geq$42M | 84.6
VGG16-GCR | 138M | 82.9
HT-LSTM (Ours) | 263K | 88.1
FDHT-LSTM (Ours) | 8,324 | 88.2
TABLE VI: Accuracy comparison with various $R$ on LSTM-UCF11 and LSTM-Youtube. Models | $R$ | Top-1 Acc.(%)
---|---|---
LSTM-UCF11 | 2 | 81.4
12 | 87.5
32 | 79.3
LSTM-Youtube | 2 | 78.9
11 | 88.2
32 | 77.4

Rank Selection. As analyzed before, for compressing large LSTM models, HT
decomposition can improve accuracy with very high compression ratio because it
can effectively reduces the structural redundancy. Therefore, for high rank
setting such as $R=32$, the accuracy is not as good as low rank setting
because considerable unnecessary redundancy still exists and it is not
friendly for training. On the other hand, for extremely low rank setting,
e.g., $R=2$, the model size is too small and thus the model capacity will be
hurt. Therefore, as shown in Table VI, in order to achieve good balance
between compression rate and model accuracy, $R=12$ and $R=11$ are set for
LSTM-UCF11 and LSTM-Youtube, respectively.

###  5.2 Performance of FDHT-LSTM Hardware

####  5.2.1 Configuration of Experiment and Design Example

A bit-accurate cycle-accurate simulator is developed to model the high-level
behavior of the proposed FDHT-LSTM architecture. Then, we build a Verilog-
based RTL model and verify the functional validation. This model is
synthesized with CMOS 28nm library. Table VII shows the detailed configuration
information of the design example based on the proposed FDHT-LSTM
architecture. Here the overall hardware consists of $16$ PEs, and each PE is
equipped with 16 16-bit multipliers and 16 24-bit accumulators, hence 256 MAC
operations are performed in each clock cycle. The working memory is composed
of $14$ banks of SRAM, and each SRAM bank has the width of 256 bits and the
depth of 2048, thereby leading to 875KB budgeted capacity that is sufficient
for storing the activation and intermediate result in most applications. Since
the working SRAM contains two copies to serve as a ping-pong buffer, the total
size of the working SRAM is $875$KB $\times 2=1.70$MB. For the weight SRAM, it
has 16-bit width with depth of 8808 to store the weight parameters of the
compressed FDHT-LSTM models. Thanks to the ultra-strong compression capability
of FDHT technique, the budgeted 17.2KB weight SRAM can support the storage of
many large-scale LSTM models.

TABLE VII: Architectural configuration for design example. System Parameter | Value
---|---
PE | 16
Quantization | 16-bit
Memory Parameter
Working Memory | Width (W) | 256-bit
Depth (M) | 2048
# of Banks (G) | 14
Total Size | 1.7MB
Weight Memory | Width | 16-bit
Depth | 8808
Size | 17.2KB
PE Parameter
Amount of Multiplier | 16
Amount of Accumulator | 16

####  5.2.2 Comparison with GPU

Though the decomposed model reduces the weight size significantly and enjoys
the parallel processing, the computation flow of FDHT-LSTM contains
complicated matrix transformation, which is not naturally supported by GPU in
an efficient way. We implement the decomposed LSTM-Youtube and LSTM-UCF11 on
NVIDIA RTX A6000, and compare the inference speed with our FDHT hardware
accelerator. The result shows that FDHT-LSTM hardware achieves 2.30$\times$
and 1.85$\times$ inference speedup than GPU for executing LSTM-Youtube and
LSTM-UCF11, respectively.

####  5.2.3 Comparison with EIE and TIE

In this subsection, we compare the proposed architecture with two prior
compressed model-oriented hardware accelerators EIE [12] and TIE [6]. Table IX
summarizes the implementation results for the three listed ASIC designs. Table
VIII shows the compression performance of FDHT-LSTM, EIE and TIE, where top-1
accuracy is the accuracy that the recognition result with the highest
probability is exactly the expected answer. It is seen that our FDHT-LSTM
outperforms EIE and TIE with respect to compression rate and accuracy
performance on both UCF11 and Youtube datasets.

TABLE VIII: Compression performance of EIE, TIE and our FDHT-LSTM on UCF11 and Youtube dataset. Top-1 accuracy is the accuracy that the recognition result with the highest probability is exactly the expected answer. | Design | Compres. Rate | Top-1 Acc.(%)
---|---|---|---
UCF11 | EIE | 117$\times$ | 82.9
TIE | 4954$\times$ | 79.6
| FDHT-LSTM
---
(Ours)
6726$\times$ | 87.5
Youtube | EIE | 111$\times$ | 84.2
TIE | 4608$\times$ | 78.8
| FDHT-LSTM
---
(Ours)
7117$\times$ | 88.2
Figure 14: Hardware performance comparison among FDHT-LSTM, TIE and EIE.

Comparison with EIE. EIE is a sparse model-oriented hardware architecture.
Consider EIE is implemented with a different technology node, as reported in
Table IX, the performance metrics of EIE are projected to 28nm technology for
fair comparison. Here the projection follows the scaling rule adopted in [12]:
linear, quadratic and constant scaling for clock frequency, silicon area and
power consumption, respectively. In addition, considering the different
budgeted arithmetic and memory resource, Figure 14 compares our proposed FDHT-
LSTM with EIE with respect to different hardware performance metrics. It is
seen that the FDHT-LSTM achieves $1.69\times$, $8.99\times$ and $6.22\times$
increase in throughput, area efficiency and energy efficiency on UCF11, and
$2.06\times$, $10.94\times$ and $7.58\times$ increase in throughput, area
efficiency and energy efficiency on Youtube respectively.

Comparison with TIE. TIE is the state-of-the-art tensor train (TT) decomposed
model-oriented hardware accelerator, which is the most related work to our
FDHT-LSTM architecture. Figure 14 shows the performance comparison between the
two designs in terms of processing throughput, area efficiency and energy
efficiency on two LSTM workloads 222Because TIE is based on the TT
decomposition that only compresses the input-to-hidden layer of the original
LSTM, the reported performance of TIE is evaluated on that single layer.. It
is seen that for executing LSTM-Youtube workload, our proposed FDHT-LSTM
architecture achieves $2.5\times$, $1.46\times$ and $2.41\times$ increase in
throughput, area efficiency and $2.41\times$ energy efficiency, respectively.
For LSTM-UCF workload, our FDHT-LSTM design also outperforms TIE with
$1.9\times$ higher throughput, $1.83\times$ higher energy efficiency and
comparable area efficiency.

TABLE IX: Comparison on hardware performance. Design | EIE [12] | TIE [6] | Our FDHT
---|---|---|---
Compression | Pruning | Tensor Train | Hierarchical
Approach | (TT) | Tucker (HT)
CMOS | 28nm | 28nm | 28nm
Technology | (projected)
Frequency | 1280 | 1000 | 1000
(MHz)
Quantization | 4-bit index | 16-bit | 16-bit
Scheme | 16-bit weight
Area (mm) | 15.7 | 1.74 | 2.96
Power (mW) | 590 | 154.8 | 160.4

####  5.2.4 Flexibility

The proposed FDHT-LSTM architecture can provide high flexibility to support
different compressed LSTM models with different topology settings. Considering
the rank $R$ is the most important parameters that directly determines the
storage requirement and computational costs of the underlying models, we study
the flexibility of the proposed hardware architecture with respect to
different $R$’s. As reported in Figure 15, FDHT-LSTM hardware demonstrates
strong flexibility.

Figure 15: Flexibility of FDHT with respect to different decomposition ranks.

##  6 Related Work

Customized hardware accelerators for deep neural networks have been
extensively studied in recent years. To accelerate the execution of
convolutional neural networks (CNNs) in a real-time and energy-efficient way,
a series of CNN hardware architecture have been proposed in [4] [24] [2]. In
particular, fully leveraging the sparsity is a very important optimization
strategy towards high-performance CNN hardware accelerators [18][1][8].

Different from computation-bounded CNNs, LSTMs are inherently memory-bounded
and require very large model size because of its component linear layer.
Therefore, the focus of efficient LSTM hardware design is to efficiently
integrate upper-level model compression algorithm to the lower-level
architecture optimization. Motivated by this design philosophy, several
compressed model-oriented LSTM hardware have been developed in
[19][11][12][6]. To be specific, EIE[12] and ESE [11] are built on pruned
model and the focus of their architectural optimization is to properly
leverage sparsity opportunities. C-LSTM [19] aims to use FFT-based compression
compression to accelerate LSTM, and hence its underlying architecture mainly
consists of optimized FFT modules. TIE is the first tensor decomposition-based
DNN hardware architecture, which is the most closed work to our design. By
leveraging the tensor train decomposition, TIE [6] can achieve competitive
hardware performance with good resource utilization.

##  7 Conclusion

In this paper, we propose algorithm and hardware co-design for FDHT-LSTM, a
ultra-compact high-performance compressed LSTM network. By fully decomposing
the entire LSTM models via hierarchical Tucker decomposition, the entire LSTM
model size can be significantly reduced. Meanwhile, an energy-efficient
customized hardware architecture with delicate design of memory access scheme
is developed. The proposed algorithm and hardware are empirically evaluated on
different video recognition datasets and LSTM workloads. The experimental
results show that our proposed FDHT-LSTM networks and the corresponding
hardware accelerator can achieve very high model performance and hardware
performance as compared to the state-of-the-art designs.

## Acknowledgments

This work was partially supported by National Science Foundation under Grant
CCF-1955909.

## References

  * [1] Alessandro Aimar, Hesham Mostafa, Enrico Calabrese, Antonio Rios-Navarro, Ricardo Tapiador-Morales, Iulia-Alexandra Lungu, Moritz B Milde, Federico Corradi, Alejandro Linares-Barranco, Shih-Chii Liu, et al.  Nullhop: A flexible convolutional neural network accelerator based on sparse representations of feature maps.  IEEE transactions on neural networks and learning systems, 30(3):644–656, 2018.
  * [2] Manoj Alwani, Han Chen, Michael Ferdman, and Peter Milder.  Fused-layer cnn accelerators.  In 2016 49th Annual IEEE/ACM International Symposium on Microarchitecture (MICRO), pages 1–12. IEEE, 2016.
  * [3] Elham Azari and Sarma Vrudhula.  An Energy-Efficient Reconfigurable LSTM Accelerator for Natural Language Processing.  In 2019 IEEE International Conference on Big Data (Big Data), pages 4450–4459, Los Angeles, CA, USA, Dec. 2019. IEEE.
  * [4] Yu-Hsin Chen, Joel Emer, and Vivienne Sze.  Eyeriss: A spatial architecture for energy-efficient dataflow for convolutional neural networks.  ACM SIGARCH Computer Architecture News, 44(3):367–379, 2016.
  * [5] Gong Cheng, Peicheng Zhou, and Junwei Han.  Duplex metric learning for image set classification.  IEEE Transactions on Image Processing, 27(1):281–292, 2017.
  * [6] Chunhua Deng, Fangxuan Sun, Xuehai Qian, Jun Lin, Zhongfeng Wang, and Bo Yuan.  TIE: energy-efficient tensor train-based inference engine for deep neural network.  In Proceedings of the 46th International Symposium on Computer Architecture, pages 264–278, Phoenix Arizona, June 2019. ACM.
  * [7] Jeffrey Donahue, Lisa Anne Hendricks, Sergio Guadarrama, Marcus Rohrbach, Subhashini Venugopalan, Kate Saenko, and Trevor Darrell.  Long-Term Recurrent Convolutional Networks for Visual Recognition and Description.  page 10.
  * [8] Ashish Gondimalla, Noah Chesnut, Mithuna Thottethodi, and TN Vijaykumar.  Sparten: A sparse tensor accelerator for convolutional neural networks.  In Proceedings of the 52nd Annual IEEE/ACM International Symposium on Microarchitecture, pages 151–165, 2019.
  * [9] Alex Graves, Navdeep Jaitly, and Abdel-rahman Mohamed.  Hybrid speech recognition with Deep Bidirectional LSTM.  In 2013 IEEE Workshop on Automatic Speech Recognition and Understanding, pages 273–278, Olomouc, Czech Republic, Dec. 2013. IEEE.
  * [10] Wolfgang Hackbusch and Stefan Kühn.  A new scheme for the tensor representation.  Journal of Fourier Analysis and Applications, 15(5):706–722, 2009\.
  * [11] Song Han, Junlong Kang, Huizi Mao, Yiming Hu, Xin Li, Yubin Li, Dongliang Xie, Hong Luo, Song Yao, Yu Wang, et al.  Ese: Efficient speech recognition engine with sparse lstm on fpga.  In Proceedings of the 2017 ACM/SIGDA International Symposium on Field-Programmable Gate Arrays, pages 75–84, 2017.
  * [12] Song Han, Xingyu Liu, Huizi Mao, Jing Pu, Ardavan Pedram, Mark A. Horowitz, and William J. Dally.  EIE: Efficient Inference Engine on Compressed Deep Neural Network.  In 2016 ACM/IEEE 43rd Annual International Symposium on Computer Architecture (ISCA), pages 243–254, Seoul, South Korea, June 2016. IEEE.
  * [13] Minyoung Kim, Sanjiv Kumar, Vladimir Pavlovic, and Henry Rowley.  Face tracking and recognition with visual constraints in real-world videos.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 1–8. IEEE, 2008.
  * [14] Yang Li, Wenming Zheng, Zhen Cui, and Tong Zhang.  Face recognition based on recurrent regression neural network.  Neurocomputing, 297:50–58, 2018.
  * [15] Bo Liu, Liping Jing, Jia Li, Jian Yu, Alex Gittens, and Michael W Mahoney.  Group collaborative representation for image set classification.  International Journal of Computer Vision, 127(2):181–206, 2019\.
  * [16] Jingen Liu, Jiebo Luo, and Mubarak Shah.  Recognizing realistic actions from videos “in the wild”.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 1996–2003. IEEE, 2009.
  * [17] Yu Pan, Jing Xu, Maolin Wang, Jinmian Ye, Fei Wang, Kun Bai, and Zenglin Xu.  Compressing recurrent neural networks with tensor ring for action recognition.  In Proceedings of the AAAI Conference on Artificial Intelligence, volume 33, pages 4683–4690, 2019.
  * [18] Angshuman Parashar, Minsoo Rhu, Anurag Mukkara, Antonio Puglielli, Rangharajan Venkatesan, Brucek Khailany, Joel Emer, Stephen W Keckler, and William J Dally.  Scnn: An accelerator for compressed-sparse convolutional neural networks.  ACM SIGARCH Computer Architecture News, 45(2):27–40, 2017.
  * [19] Shuo Wang, Zhe Li, Caiwen Ding, Bo Yuan, Qinru Qiu, Yanzhi Wang, and Yun Liang.  C-lstm: Enabling efficient lstm using structured compression techniques on fpgas.  In Proceedings of the 2018 ACM/SIGDA International Symposium on Field-Programmable Gate Arrays, pages 11–20, 2018.
  * [20] Yinchong Yang, Denis Krompass, and Volker Tresp.  Tensor-train recurrent neural networks for video classification.  In International Conference on Machine Learning, pages 3891–3900. JMLR. org, 2017.
  * [21] Jinmian Ye, Linnan Wang, Guangxi Li, Di Chen, Shandian Zhe, Xinqi Chu, and Zenglin Xu.  Learning compact recurrent neural networks with block-term tensor decomposition.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 9378–9387, 2018.
  * [22] Miao Yin, Siyu Liao, Xiao-Yang Liu, Xiaodong Wang, and Bo Yuan.  Towards extremely compact rnns for video recognition with fully decomposed hierarchical tucker structure.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 12085–12094, 2021.
  * [23] Rose Yu, Stephan Zheng, Anima Anandkumar, and Yisong Yue.  Long-term forecasting using higher order tensor rnns.  arXiv preprint arXiv:1711.00073, 2017.
  * [24] Chen Zhang, Peng Li, Guangyu Sun, Yijin Guan, Bingjun Xiao, and Jason Cong.  Optimizing fpga-based accelerator design for deep convolutional neural networks.  In Proceedings of the 2015 ACM/SIGDA international symposium on field-programmable gate arrays, pages 161–170, 2015.

|  Yu Gong received the master’s degree in Shanghai Jiao Tong University, and
bachelor’s degree in Wuhan University of Technology. Currently, he is a Ph.D.
student in Electrical and Computer Engineering at Rutgers University. His
research interests include high-performance architecture for AI and hardware-
software co-design.
---|---
|  Miao Yin is a Ph.D. student in Electrical and Computer Engineering at
Rutgers, The State University of New Jersey. His research is focused on
algorithm-hardware co-designed energy-efficient AI system (e.g., on-device
image and action recognition) based on higher-order tensor decomposition and
advanced optimization. He serves as a PC member of the IEEE/CVF Conference on
Computer Vision and Pattern Recognition.
---|---
|  Lingyi Huang received the master’s degree in electrical engineering from
New York University and the bachelor’s degree in electrical engineering from
Xidian University. He is currently pursuing the PhD degree in computer
engineering at Rutgers University.
---|---
|  Chunhua Deng is a senior staff design engineer at ScaleFlux Inc. He
received his Ph.D. degree from the Department of Electrical and Computer
Engineering at Rutgers University. He received his bachelor’s degree and
master’s degree from China University of Petroleum, Beijing Institute of
Technology, respectively. His research interests include machine learning,
computer architecture, and VLSI design.
---|---
|  Yang Sui is a Ph.D. student in Electrical and Computer Engineering at
Rutgers, The State University of New Jersey. His research is focusing on
algorithm-hardware co-designed efficient AI with model compression based on
pruning, low-rank decomposition, and advanced algorithms.
---|---
|  Bo Yuan received the bachelor’s and master’s degrees from Nanjing
University, Nanjing, China, in 2007 and 2010, respectively, and the PhD degree
from the Department of Electrical and Computer Engineering, University of
Minnesota, Twin Cities, Minnesota, in 2015. His research interests include
algorithm and hardware co-design and implementation for machine learning and
signal processing systems, error-resilient low-cost computing techniques for
embedded and IoT systems and machine learning for domain-specific
applications. He is the recipient of Global Research Competition Finalist
Award in Broadcom Corporation and doctoral dissertation fellowship with the
University of Minnesota. He serves as technical committee track chair and
technical committee member for several IEEE/ACM conferences. He is the
associated editor of the Springer Journal of Signal Processing System.
---|---


# A physics-informed search for metric solutions to Ricci flow, their
embeddings, and visualisation

Aarjav Jain
Dept. of Physics
Cambridge University
<EMAIL_ADDRESS>
&Challenger Mishra
Dept. of Computer Science
Cambridge University
cm2099@@cam.ac.uk
&Pietro Liò
Dept. of Computer Science
Cambridge University
<EMAIL_ADDRESS>

###### Abstract

Neural networks with PDEs embedded in their loss functions (physics-informed
neural networks) are employed as a function approximators to find solutions to
the Ricci flow (a curvature based evolution) of Riemannian metrics. A general
method is developed and applied to the real torus. The validity of the
solution is verified by comparing the time evolution of scalar curvature with
that found using a standard PDE solver, which decreases to a constant value of
0 on the whole manifold. We also consider certain solitonic solutions to the
Ricci flow equation in two real dimensions. We create visualisations of the
flow by utilising an embedding into $\mathbb{R}^{3}$. Snapshots of highly
accurate numerical evolution of the toroidal metric over time are reported. We
provide guidelines on applications of this methodology to the problem of
determining Ricci flat Calabi–Yau metrics in the context of String theory, a
long standing problem in complex geometry.

##  1 Introduction

Ricci flow is a powerful tool for geometric analysis [1]. Formulated by
Hamilton in the 80’s, Grigori Perelman employed Ricci flow to resolve the
Poincaré conjecture in 2002, proving that every simply connected, closed
manifold in three dimensions is homeomorphic to the 3-sphere [2]. Intuitively,
Ricci flow is a geometric flow of a metric on a smooth Riemannian manifold,
which uniformizes curvature of the metric. Given a smooth Riemannian manifold
$\mathcal{M}$, Ricci flow assigns a metric _g_ to each _time_ $t$, dictated by
the Ricci tensor of the metric, Ric(g) [3]. Mathematically put,

$\partial_{t}g_{ab}(x,t)=-2Ric_{ab}(g(x,t)),\leavevmode\nobreak\
\text{where\leavevmode\nobreak\ $a$,
$b$}\in\\{1,2,\dots,dim({\mathcal{M}})\\}.$ (1)

Ricci flow has deep connections to renormalization group flow, the latter
appearing in the context of beta functions in quantum field theories, which
formulate the dependency of coupling constants on energy scales [4, 5]. The
flow parameter $t$ in Equation 1 is reminiscent of such scales. Following
decades of further advances, Ricci flow has found connections with training in
neural networks, as well as in community learning [6] and studying robustness
under gradient descent perturbations [7].

Our interest in Ricci flow concerns complex geometries which appear as extra
dimensions in certain superstring theories. A class of geometries of
particular interest is Calabi–Yau manifolds. One problem of particular
interest is related to the Calabi conjecture, which postulated the existence
of metrics over such geometries that have a vanishing Ricci tensor, that is,
‘Ricci flat’. This amounts to finding the fixed point of the PDE in Equation
1\- an analytic form of which remains elusive till date. Further work in the
interim resulted in numerical treatments leading to approximation techniques
[8] including the celebrated Donaldson’s algorithm [9]. Recently, machine
learning was demonstrated to be a powerful tool in the context of
approximating Ricci flat Calabi–Yau metrics for complex geometries with
special holonomy [10, 11, 12]. These approaches rely on neural network
approximations of the metric of interest. This is reminiscent of physics
inspired neural networks (PINNS)[13]; albeit the only data in the machine
learning pipeline involves points sampled from such geometries. Hinted by
[10], in this work we propose to formulate Ricci flow on real geometries using
PINNs. We focus on the real torus, and consider cigar solitonic solutions to
the Ricci flow equation 1.

##  2 Using neural networks to learn Ricci flow

The flow described in 1, is the ‘unnormalised’ Ricci flow. As such the volume
of the given Riemannian manifold can evolve with time. The action of the flow
is to uniformize curvature over the surface, in a way analogous to the heat
equation spreading heat evenly over space. The action of Ricci flow on the
round metric of a sphere can be shown to shrink the sphere to a point [3]. We
utilise the Nash–Kuiper theorem to embed the geometries in $\mathbb{R}^{3}$.
We propose to use Physics-informed neural networks (PINNs) to solve the Ricci
flow PDE for two-dimensional real manifolds. PINNs incorporate the PDEs, and
the key symmetry/initial conditions of the system into the loss function,
meaning that once trained, the neural network will output a solution to the
PDE ([14]). The metric $g{(x,t)}$, on the manifold $\mathcal{M}$, is
approximated by the neural network $g_{\theta}(x,t)$, where $\theta$ denotes
the parameters of the neural network. The PDE residual is calculated from the
Ricci flow equation 2 using the neural network approximation. From this, we
define a residual loss function, $\mathcal{L}_{res}{(\theta)}$ in Equation 2.
Similarly, residual losses associated with the required initial and symmetry
conditions are constructed. These are $\mathcal{L}_{init}{(\theta)}$ and
$\mathcal{L}_{sym}{(\theta)}$ respectively, and they are combined to make the
PINN loss function in Equation 3. The $\mathcal{L}_{sym}{(\theta)}$ loss can
include terms which mandate a desired symmetry. If we know the initial metric
has rotational symmetry, we could enforce this by inclusing a rotational
symmetry loss function. More on this in Section 3.

$\displaystyle\mathcal{L}_{res}{(\theta)}:=\frac{1}{N}\sum^{N}_{i=1}\left|r_{\theta}(t_{i},x_{i})\right|^{2},\text{\leavevmode\nobreak\
with\leavevmode\nobreak\ \leavevmode\nobreak\
}r_{\theta}(t,x):=\frac{\partial}{\partial
t}g_{\theta}(t,x)+2Ric(g_{\theta}(t,x))\text{\leavevmode\nobreak\
,\leavevmode\nobreak\ }x_{i}\in\mathcal{M}.$ (2)
$\mathcal{L}{(\theta)}=\mathcal{L}_{res}{(\theta)}+\lambda_{0}\
\mathcal{L}_{\text{init}}{(\theta)}+\lambda_{1}\
\mathcal{L}_{\text{sym}}{(\theta)}.$ (3)

##  3 Experiments

We apply the PINN approach to Ricci flow of a metric on the real 2-dimensional
torus. We also consider cigar solitonic solutions of Ricci flow as an
additional study, since this has analytic solutions to Ricci flow, hence
allows the benchmarking of PINN generated solutions.

Real Torus:  The torus is spanned by toroidal and poloidal angles, denoted by
$u$ and $v$ respectively, with $u\in[0,2\pi)$ and $v\in[0,2\pi)$. An
illustration of this parameterisation is in the Figure 5 in the Appendix A. A
natural metric tensor on the real torus is,

$\displaystyle g(t,u,v)=\begin{pmatrix}(c+a\cos v)^{2}&0\\\
0&a^{2}\end{pmatrix},\text{\leavevmode\nobreak\ \leavevmode\nobreak\ with
symmetry conditions, }$ (4) $\displaystyle
g(t,u,v)=g(t,u+\delta,v)=g(t,u,2\pi-v),\text{ where }\delta\in[0,2\pi),$ (5)

which follow from the symmetries of the torus. Training is done using the Adam
optimiser [15] on a fully connected feed-forward network with three hidden
layers of 16, 32 and 16 units, and Softplus activation. We use 1000 randomly
generated training points in the domain. Our implementation is in JAX which
supports efficient gradient computations, through auto differentiation. This
is particularly important in light of the fact that the Ricci tensor,
$Ric(g)$, is dependent on the derivatives of the metric function. In Figures 1
and 2, we show the evolution of the scalar curvature (which is the trace of
the $Ric(g)$) with the flow parameter. We compare the PINN generated solution
with numerical simulations using a standard PDE solver in Mathematica. We
choose $c$ = 2 and $a$ = 1. The solutions bear significant resemblance to each
other. The initial Ricci scalar curve matches the exact result which can be
found by taking derivatives of the initial torus metric shown in Equation 4\-
a useful sanity check. The PINN and exact results have an MSE of $\sim
1.8\times 10^{-3}$. The computation took around 4 hours on an i5 laptop. An
additional visualisation is provided in Figure 2, where the torus base
geometry is kept the same, but the decay of curvature with time under Ricci
flow is depicted with colour representing scalar curvature. Curvature decays
to a zero with time, as expected.

(a)

(b)

Figure 1: Evolution of the scalar curvature under Ricci flow. A Mathematica
PDE solver generated the left figure; the right corresponds to our PINN
solution. It captures the initial conditions, and shows accurate long time
behaviour of vanishing scalar curvature, the fixed point of Ricci flow.

Figure 2: Colour plot of Ricci scalar of a torus evolving under Ricci flow.
Here, unlike the embedding case, the base geometry is kept the same,
showcasing the decay of curvature with Ricci flow.

Figure 3: Left is a plot of the initial torus metric. Its geodesics are
plotted in colours. Right is a plot of the torus Ricci flow ($t=2$). Note the
scale in the $z$ direction. It does not have the required rotational symmetry,
but has been severely flattened in the $z$ direction, as expected.

(a)

(b)

Figure 4: The figure on the left is the $g_{11}=g_{22}$ component of the Cigar
soliton solution generated by PINN under Ricci flow. This compares favourably
with the analytic solution in Equation 6. On the right we show the MSE. It is
$\mathcal{O}(10^{-2})$ at small times, and grows under Ricci flow.

Cigar Soliton:  An example of a two dimensional metric with an analytic
solution to the Ricci flow is the cigar soliton [16]. It has metric
$g(\textbf{x})=({1+x^{2}+y^{2}})^{-1}\mathbb{I}_{2}$, with its Ricci flow
solution given by:

$g(t,\textbf{x})=({e^{4t}+x^{2}+y^{2}})^{-1}\mathbb{I}_{2},\text{\leavevmode\nobreak\
where,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
}(x,y)\in\mathbb{R}^{2}.$ (6)

The PINN solution approximates the analytic metric for all values of the flow
parameter, $t$. A three dimensional plot of the $g_{11}=g_{22}$ component as a
function of $x$ and $y$ is shown on the left of Figure 4. The behaviour of the
PINN generated solution bears close resemblance to the analytic solution,
although over time, the MSE gradually grows. We expect high accuracies at
lower times since the residuals from the metric at $t=0$ is minimised by
choice of loss function (Equation 3). The minima occurs since, as time
increases, the metric values themselves decrease, and this is overcome at
higher time.

##  4 Metric plotting

We aim to plot the surface corresponding to the learned neural network
solution to Ricci flow, $g_{\theta}(t,u,v)$. That is, we wish to create an
embedding of the manifold into $\mathbb{R}^{3}$, and plot this for a number of
randomly generated points to give a three dimensional plot of the manifold for
multiple time durations under Ricci flow. A smooth embedding of the geometries
we have considered into Euclidean space exists ala the Nash-Kuiper embedding
theorem. We denote this by $f:M\rightarrow\mathbb{R}^{n}$. The embedding is
given by a system of first-order PDEs:

$g_{ij}(x)=\sum_{\alpha=1}^{n}{\frac{\partial f^{\alpha}}{\partial
x^{i}}}{\frac{\partial f^{\alpha}}{\partial x^{j}}},\alpha\in 1,2,\dots,n$ (7)

where index $\alpha$ runs over the number of dimensions for the target
embedding (3 in this case). Embedding into $\mathbb{R}^{3}$ leaves a system of
the three PDEs (An explicit unpacking is detailed in Appendix A), the solution
of which give $f^{i}(u,v)$, where $i=dim(\mathcal{M})$. This embedding PDE is
solved using PINNs. The architecture, epochs and training splits are identical
to Section 3. The training data consists of 1000 data points with random
angles on the torus at a fixed time point, with the output being the
corresponding points in $\mathbb{R}^{3}$ (embedding). We proceed by sampling
random points on the base metric $(u_{i},v_{i})$ and applying the embedding
functions $f^{i}$ to each of these, as

$(x_{i},y_{i},z_{i})=(f^{x}(u_{i},v_{i}),f^{y}(u_{i},v_{i}),f^{z}(u_{i},v_{i})).$
(8)

These points in $\mathbb{R}^{3}$ can be plotted to give a visualisation of the
manifold at different times during Ricci flow. Ricci flow is diffeomorphism
invariant, so the symmetries of the initial metric should be preserved. Hence,
symmetries in Equation 5 are included in the loss function through
$\mathcal{L}_{\text{sym}}$ in Equation 2. By embedding 1000 random points on
the torus into $\mathbb{R}^{3}$ and plotting them we get plots for initial
torus, and torus after 2 seconds of Ricci flow, seen in Figure 3. The geodesic
equation is solved for trivial initial directions using this embedding, and
the curves are plotted on the surface of the initial torus, giving the
expected equatorial paths. Training took 8 hours on an i5 laptop. After two
seconds, the manifold is much flatter, which is expected (notice the scale in
$z$-direction). In the figure on the right, the expected rotational symmetry
is not seen. It was found that this plot depended heavily on the
regularisation terms applied to various loss terms (see Equation 2). In this
experiment, fixed scalars were used for these terms, but if an adaptive scheme
is used (such as [17]), then a term corresponding to rotational symmetry could
be enhanced when necessary during training.

##  5 Discussion and Limitations

We applied a PINN technique to solve the Ricci Flow PDE for real, two-
dimensional geometries. The solution for a torus matched well to that found by
a PDE solver, and had an MSE of $\mathcal{O}(10^{-3})$ compared with the exact
solution for the toroidal initial metric. The cigar soliton geometry matched
the analytic solution, with an MSE of $\mathcal{O}(10^{-2})$. The learned PINN
metric for torus Ricci flow was used to solve PDEs derived from the
Nash–Kuiper Embedding theorem using PINNs. This allows visualisations of
manifolds under Ricci Flow, where the geodesics of the metric are considered
as spanning the base geometry. These showed elements of correctness; a
toroidal figure was seen before Ricci Flow, and a flatter manifold after Ricci
Flow. Full rotational symmetry was not perfectly replicated in either case,
which are potentially due to compounding numerical errors during the course of
Ricci flow, although the growth shows signs of stabilising. Instead of
demanding rotational symmetry from the loss function, one may embed it
directly with an appropriate coordinate system. This was avoided here, since
choosing one can become difficult when working with larger dimensional and/or
complex manifolds. Although deliberate in our choice of architecture, the
approach would benefit from the use of equivariant architectures. Further, the
magnitudes of various loss terms being important to the final embedding, an
adaptive technique such as [17] could be employed. The PINN approach allowed
simple differentiation of the neural network solution, to find relevant
quantities of Ricci tensor and curvature, as well as being beneficial in the
context of Nash–Kuiper embedding. The results presented show the promise of
using PINN techniques in Ricci flow of higher dimensional manifolds and
complex manifolds, and could aid in their visualisation. The tool we have
created is completely general, and allows adaptation to the setting of real
manifolds by simply changing the initial metric loss function.

Embeddings of complex analogue of the real torus, K3 surfaces, and Calabi–Yau
manifolds are possible, when technical aspects related to sampling and
constructing meaningful loss functions [11, 12, 10] are considered. See A for
further details. Solutions to Ricci flow have substantial applications in
superstring theory wherein extra dimensional geometries play a crucial role in
determining four-dimensional physics, enabling predictions of particle masses
and coupling constants from string theory for the first time, which rely on an
accurate computation of the Ricci flat metric.

##  6 Broader Impact

Ricci flow is mathematically rich and a powerful tool for geometric analysis
We develop and demonstrate a framework for learning solutions to Riemannian
metrics, and provide special cases for real 2$D$ geometries. This framework
can be extended to complex geometries with increasing dimension, such as the
complex torus, K3 surfaces, and Calabi–Yau manifolds, which have applications
in superstring theory: enabling prediction of particle masses. We also
presenet an embedding framework of geometries through PINNS and this idea can
aid visualisation of higher dimensional manifolds.

## References

  * [1] Bennett Chow, Peng Lu and Lei Ni  “Hamilton’s Ricci flow”  American Mathematical Soc., 2006
  * [2] Grisha Perelman  “The entropy formula for the Ricci flow and its geometric applications”  In _arXiv preprint math/0211159_ , 2002
  * [3] Nick Sheridan, Craig Hodgson and J. Rubinstein  “Hamilton’s Ricci Flow”, 2006
  * [4] Ioannis Bakas  “Renormalization group equations and geometric flows”  In _arXiv preprint hep-th/0702034_ , 2007
  * [5] Mauro Carfora  “Renormalization group and the Ricci flow”  In _Milan journal of mathematics_ 78.1 Springer, 2010, pp. 319–353
  * [6] Samuel Glass, Simeon Spasov and Pietro Liò  “RicciNets: Curvature-guided Pruning of High-performance Neural Networks Using Ricci Flow”  In _arXiv preprint arXiv:2007.04216_ , 2020
  * [7] Jun Chen, Tianxin Huang, Wenzhou Chen and Yong Liu  “Thoughts on the Consistency between Ricci Flow and Neural Network Behavior”  In _arXiv preprint arXiv:2111.08410_ , 2021
  * [8] Matthew Headrick and Toby Wiseman  “Numerical Ricci-flat metrics on K3”  In _Classical and Quantum Gravity_ 22.23 IOP Publishing, 2005, pp. 4931
  * [9] S.. Donaldson  “Some numerical results in complex differential geometry”  arXiv, 2005  DOI: 10.48550/ARXIV.MATH/0512625
  * [10] Vishnu Jejjala, Damian Kaloni Mayorga Pena and Challenger Mishra  “Neural network approximations for Calabi-Yau metrics”  In _Journal of High Energy Physics_ 2022.8 Springer, 2022, pp. 1–37
  * [11] Michael Douglas, Subramanian Lakshminarasimhan and Yidi Qi  “Numerical Calabi-Yau metrics from holomorphic networks”  In _Mathematical and Scientific Machine Learning_ , 2022, pp. 223–252  PMLR
  * [12] Lara B Anderson et al.  “Moduli-dependent Calabi-Yau and SU (3)-structure metrics from Machine Learning”  In _Journal of High Energy Physics_ 2021.5 Springer, 2021, pp. 1–45
  * [13] Maziar Raissi, Paris Perdikaris and George Em Karniadakis  “Physics informed deep learning (part i): Data-driven solutions of nonlinear partial differential equations”  In _arXiv preprint arXiv:1711.10561_ , 2017
  * [14] Salvatore Cuomo et al.  “Scientific Machine Learning through Physics-Informed Neural Networks: Where we are and What’s next”  arXiv, 2022  DOI: 10.48550/ARXIV.2201.05624
  * [15] Diederik Kingma and Jimmy Ba  “Adam: A Method for Stochastic Optimization”  In _International Conference on Learning Representations_ , 2014
  * [16] Peter Topping  “Lectures on the Ricci Flow”, London Mathematical Society Lecture Note Series  Cambridge University Press, 2006  DOI: 10.1017/CBO9780511721465
  * [17] Zixue Xiang, Wei Peng, Xu Liu and Wen Yao  “Self-adaptive loss balanced Physics-informed neural networks”  In _Neurocomputing_ 496, 2022, pp. 11–34  DOI: https://doi.org/10.1016/j.neucom.2022.05.015
  * [18] Philip Candelas, Anders Martin Dale, CA Lütken and Rolf Schimmrigk  “Complete intersection calabi-yau manifolds”  In _Nuclear Physics B_ 298.3 Elsevier, 1988, pp. 493–525
  * [19] Tristan Hubsch  “Calabi-Yau manifolds: A Bestiary for physicists”  World scientific, 1992
  * [20] P. Candelas, C.. Lutken and R. Schimmrigk  “COMPLETE INTERSECTION CALABI-YAU MANIFOLDS. 2. THREE GENERATION MANIFOLDS”  In _Nucl. Phys. B_ 306, 1988, pp. 113  DOI: 10.1016/0550-3213(88)90173-3
  * [21] Philip Candelas and Rhys Davies  “New Calabi-Yau Manifolds with Small Hodge Numbers”  In _Fortsch.Phys._ 58, 2010, pp. 383–466  DOI: 10.1002/prop.200900105
  * [22] Philip Candelas and Andrei Constantin  “Completing the Web of $Z_{3}$ \- Quotients of Complete Intersection Calabi-Yau Manifolds”  In _Fortsch.Phys._ 60, 2012, pp. 345–369  DOI: 10.1002/prop.201200044
  * [23] Volker Braun  “On Free Quotients of Complete Intersection Calabi-Yau Manifolds”  In _JHEP_ 1104, 2011, pp. 005  DOI: 10.1007/JHEP04(2011)005
  * [24] Andre Lukas and Challenger Mishra  “Discrete symmetries of complete intersection Calabi–Yau manifolds”  In _Communications in Mathematical Physics_ 379.3 Springer, 2020, pp. 847–865
  * [25] Challenger Mishra  “Calabi-Yau manifolds, discrete symmetries and string theory”, 2017
  * [26] Philip Candelas and Challenger Mishra  “Highly symmetric quintic quotients”  In _Fortschritte der Physik_ 66.4 Wiley Online Library, 2018, pp. 1800017

##  Appendix A Appendix

In this appendix we provide some further details on the parameterisation of
the real torus, and an explicit form of the partial differential equations
that appear in the Nash-Kupier Embedding theorem.

###  A.1 Embeddings and Ricci flow

Figure 5: A diagram depicting the the toroidal ($u$) direction, represented by
the blue arrow, and the poloidal ($v$) direction, represented by the red arrow

#### Torus curvature decay

During Ricci flow, the scalar curvature $R$ of the torus is expected to reduce
with time. Figure 2 shows stills of an animation produced to showcase this.
The animation can be found in the GitHub link.

#### Nash-Kupier Embedding equation

$g_{ij}(x)=\sum_{\alpha=1}^{n}{\frac{\partial f^{\alpha}}{\partial
x^{i}}}{\frac{\partial f^{\alpha}}{\partial x^{j}}},\alpha\in 1,2,\dots,n$ (9)

###  A.2 Complete Intersection Calabi–Yau manifolds

Calabi–Yau manifolds are complex, Kähler manifolds with vanishing first Chern
class, or equivalently, with holonomy group $SU(n)$, where $n$ is the complex
dimension of the manifold. A large class of Calabi–Yau manifolds, is the
Complete Intersection Calabi–Yau manifolds (CICYs) [18]. Compiled in the late
80’s, this is the first such dataset in algebraic-geometry. The ambient space
for these geometries is a product of complex projective spaces (as opposed to
$\mathbb{R}^{n}$). The embedding is given by the zero locus of intersecting
polynomials with degrees that satisfy certain constraints that ensure the
Calabi–Yau property [19]. The class of Complete Intersection Calabi–Yau
manifolds (CICYs) in three dimensions are a family of 7890 geometries. A CICY
$\mathcal{M}$ is described by the vanishing locus of a number of polynomials
with co-ordinates in an ambient space $\cal{A}$ which is a product of complex
projective spaces of different dimensions. That is,

$\mathcal{M}\subset{\cal
A}=\mathbb{P}^{n_{1}}\times\cdots\times\mathbb{P}^{n_{m}}\;.$

The deformation class of $\mathcal{M}$ is given by a configuration matrix
which contains the multi-degrees the defining polynomials. These geometries
have an equivalent bipartite graph representation in which the polynomials and
components of the projective spaces are represented as distinct types of
nodes. The multi-degrees are captured by the number of connections between a
given space and a polynomial. A couple of examples follow.

The bipartite graph on the left, refers to a Calabi–Yau geometry, known as
tetraquadric. The one to the right is a split of this [19], which is a
geometric transition between the two manifolds. It is apparent from these
graphical representations that these geometries can admit discrete symmetries
[20, 21, 22, 23, 24, 25, 26, 25] and as such are amenable to a treatment
developed in this paper using equivariant architectures.