# QuadraNet: Improving High-Order Neural Interaction Efficiency
with Hardware-Aware Quadratic Neural Networks

Chenhui Xu*, Fuxun Yu*, Zirui Xu*, Chenchen Liu$\dagger$, Jinjun
Xiong$\ddagger$, Xiang Chen$*$,# *George Mason University, $\dagger$University
of Maryland, $\ddagger$University at Buffalo, #Peking University
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Recent progress in computer vision-oriented neural network designs is mostly
driven by capturing high-order neural interactions among inputs and features.
And there emerged a variety of approaches to accomplish this, such as
Transformers and its variants. However, these interactions generate a large
amount of intermediate state and/or strong data dependency, leading to
considerable memory consumption and computing cost, and therefore compromising
the overall runtime performance. To address this challenge, we rethink the
high-order interactive neural network design with a quadratic computing
approach. Specifically, we propose QuadraNet — a comprehensive model design
methodology from neuron reconstruction to structural block and eventually to
the overall neural network implementation. Leveraging quadratic neurons’
intrinsic high-order advantages and dedicated computation optimization
schemes, QuadraNet could effectively achieve optimal cognition and computation
performance. Incorporating state-of-the-art hardware-aware neural architecture
search and system integration techniques, QuadraNet could also be well
generalized in different hardware constraint settings and deployment
scenarios. The experiment shows that QuadraNet achieves up to 1.5$\times$
throughput, 30% less memory footprint, and similar cognition performance,
compared with the state-of-the-art high-order approaches.

##  I Introduction

During the past decade, backbone neural networks for computer vision tasks
have evolved from convolution-based information extraction [1, 2] into
Transformer-based semantic representation that well leverages self-attention
mechanisms for information superposition [3, 4]. However, such an advanced
mechanism brings huge computation costs. Specifically, self-attention
generation involves complex neural information interactions as shown in Fig. 1
(a), resulting in more computing resource consumption and memory overheads.
Therefore, even with supreme cognition performance, an extremely lightweight
MobileViT model [5] for edge devices still suffers a 10$\times$ latency slow-
down than a CNN-based MobileNet_V2 [6], among which more than 60% slow-down is
caused by the resource-consuming self-attention generation.

As self-attention generation emerges as a new computing bottleneck for
Transformer-like neural networks, an open question is how to better understand
its in-depth working mechanism. From the algorithmic perspective, self-
attention generations are achieved by high-order neural interactions of
iterative point products across multiple transformed feature dimensions,
generating input-adaptive weight parameters. Although such high-order neural
interactions may best be captured at the cost of quadratic time and space
complexity [7], a large number of sparse intermediate states are introduced as
shown in Fig. 1 (b). Compared to regular feature convolution operations with a
linear complexity to the input size, these intermediate states consume a
quadratic feature space, thus resulting the aforementioned high computation
cost. From the computing perspective, to avoid self-attention costs, some
recent works proposed other means to achieve similar high-order neural
interactions. For example, [8] bypassed softmax operations to approximate
self-attention as associative multiplication to avoid the quadratic
intermediate states. [9] also devolved high-order production with gated
recursive convolutions for lower computation complexity. However, such
approaches suffer from a number of problems. For example, many self-attention
approximations fell into a sub-optimal feature space, compromising the
cognition performance [8]; or some achieved comparable cognition performance
but introduced more intricate data dependency, limiting not only the inference
throughput but the model design generalization [9].

Figure 1:  High-Order Neural Interactions in Self-Attention Generations

Different from these state of arts [9, 8], in this paper, we propose to
rethink the Transformer-like neural network constructions’ in capturing the
high-order neural interactions. Based on that understanding, we then propose a
new neuron-level design coupled with neural architecture designs to address
the high-order issues while retaining high accuracy, low computation, and
great generalization capacities. Inspired by [10], which builds neural
networks with polynomials without activation functions, we found that high-
order neural interactions could be directly embedded into high-order neuron
designs. In particular, such neuron implementations could be achieved by
quadratic neurons [11, 12]. Quadratic neurons possess a strong ability for
information self-reinforcement across multiple feature dimensions with much
simpler and more efficient computation patterns, therefore eliminating
intermediate states from conventional self-attention generations. Thus, we
take quadratic neurons as the key to the revolution of fundamental
Transformer-like model designs.

Following such an approach, in this paper, we establish a novel quadratic
neural network architecture — QuadraNet, which masters high-order neural
interactions and adapts to different hardware computing constraints. To take
inspiration from quadratic neurons to a successful comprehensive QuadraNet
architecture realization, we need to tackle a series of challenges, including
designing a more efficient quadratic neuron and leveraging emerging compiling
techniques to accelerate quadratic convolution operations; building basic
model structural blocks by incorporating quadratic convolutions into modern
network models; and incorporating hardware awareness into the design of high-
order interactive quadratic networks through a neural architecture search
(NAS) approach.

The main contributions of this paper are as follows.

  * •

We demystify the high-order neural interactions in Transformer-like neural
networks through comprehensive analyses of algorithms and computation
patterns, and explore the bottleneck of existing optimization approaches.

  * •

We break through conventional Transformer-like model designs with a full-stack
innovation from quadratic neurons to structural blocks and eventually to
comprehensive artificially designed and searched neural networks.

  * •

We further generalize the quadratic neural networks design methodology to
various model architectures and heterologous hardware computing constraints.

To the best of our knowledge, QuadraNet is the very first quadratic neural
network design methodology that achieves competitive cognition performance
with state-of-the-art Transformer-like models while achieving outstanding
computation efficiency. Specifically, the experiment shows that the proposed
QuadraNet performs exceptional high-order neural interactions in terms of
comparative self-attention representation capacity and classification
accuracy. Meanwhile, QuadraNet can provide 1.5$\times$ throughput without any
accuracy reduction and 30% memory occupation optimization compared with
representative Transformer-like models [4].

##  II Motivation Analysis

###  II-A High-Order Interactive Neural Networks

Existing Approaches: To capture high-order neural interactions, Transformer-
like neural network models, such as ViT [3] and its variants [4], adopt the
self-attention mechanism. In conducting model cognition, self-attention
mechanism performs a spatial mixing computation $y_{SA}^{(i,c)}$ of the
$c^{th}$ output channel at the $i^{th}$ input feature dimension:

$y_{Self-
Attention}^{(i,c)}=\sum_{j\in\Omega_{i}}\sum_{c^{\prime}=1}^{C}A(x)_{i\rightarrow
j}T_{V}^{(c^{\prime},c)}x^{(j,c^{\prime})},$ (1)

where $\Omega_{i}$ is the attention receptive field involving interactive
neurons; $A(x)$ is the input-adaptive weight matrix calculated by a softmax
dot-production of “query” ($Q$) and “key” ($K$) transformations of the input
$x$; and the $i$ $\rightarrow$ $j$ process denotes the relative position of
the $i^{th}$ and $j^{th}$ feature dimension. Meanwhile, self-attention
achieves channel-wise information transformation through its “value” ($V$)
transformation of $T_{V}$ on $x^{(j,c^{\prime})}$, and $(c^{\prime},c)$
denotes the information transformation as manifested from the $c^{\prime th}$
input channel to the $c^{th}$ output channel. We argue that such a
multiplication of $A(x)$ and $x$ forms high-order neural interactions, which
is the key to such neural model’s high cognition performance.

Different from the self-attention shown in Eq. (1), HorNet [9] performs high-
order neural interactions with a gated recursive convolution — grConv, which
can be formulated as:

$y_{g^{r}Conv}^{(i,c)}=p_{r}^{(i,c)}=\sum_{j\in\Omega_{i}}\sum_{c^{\prime}=1}^{C}g(p_{r-1})w_{r-1,i\rightarrow
j}^{c}T_{\phi_{in}}^{(c^{\prime},c)}x^{(j,c^{\prime})}.$ (2)

Instead of generating high-order interactions directly from the input $x$,
grConv uses a channel increasing recursive formulation $g(p_{r-1})$ to perform
a high-order iteration of feature information, and thus getting an input-
adaptive weight matrix. $T_{\phi_{in}}^{(c^{\prime},c)}$ is a linear dimension
switch input projection.

Higher-Order Neural Interaction: By summarizing the above works, we
conclusively find that both approaches adopt an input-dependent weight matrix
and a channel-wise information transformation matrix. It can be formulated as
a generic high-order neural interaction mechanism:

$y_{high-
order}^{(i,c)}=\sum_{j\in\Omega_{i}}\sum_{c^{\prime}=1}^{C}g(x)_{ij}T^{(c^{\prime},c)}x^{(j,c^{\prime})},$
(3)

where $\Omega_{i}$ denotes the mechanism’s receptive field; $g(x)$ denotes an
interactive weight matrix that contains the input’s information; and
$T^{(c^{\prime},c)}$ denotes a channel-wise information transformation. The
multiplication of $x$ and input-adaptive $g(x)$ forms the high-order neural
interaction of data itself. Such an interaction leverages the superposition of
activation information between different feature dimensions.

###  II-B High-Order Computation Problems

The computation of existing high-order interaction approaches (esp., self-
attention) face multiple problems regarding memory consumption and computing
cost.

Memory Consumption with Intermediate States: Conventional Transformer-like
models usually suffer from excessive memory consumption. For example, in the
training stage, all the data in intermediate states must be saved for
backpropagations; in the inference stage, frequent and continuous generation
and release of intermediate states during forward propagation also cast
inconsistency and interval of memory operation and thus heavy IO consumption
in runtime. Such consumption significantly restricts the neural model
architecture’s design flexibility and optimization feasibility.

As aforementioned, those intermediate states are introduced by self-attention
generation. Specifically, the transformation process of $Q$, $K$, and $V$
parameters in Transformer-like models generates 3$HWC$ intermediate states,
and the $softmax$($QK^{T}$)$V$ takes intermediate states as many as
($HW$)2+($HW$)${}^{2}C$, where $H$, $W$ and $C$ denote the height, width, and
channels of the input per layer, respectively. Compared with widely-used
depthwise separable convolution, which generates only $HWC$ intermediate
states, self-attention needs $\sim$($HW$+3)$\times$ for computation. Even if
the attention space is restricted to a local receptive field as in [4], the
self-attention computation still yields $\sim$($M^{2}$+3)$\times$ more
intermediate states, where $M$ denotes the field size (usually $\geq$7 [4]).

Computing Cost with Data Dependency: In addition to memory consumption, the
complex computing flow in self-attention also introduces huge costs. For
example, the computation of grConv in HorNet makes use of recursive
convolutions to achieve higher-order neural interactions. Although such
ingenious design limits the FLOPs to an acceptable range, the recursive logic
brings considerable data dependency issues. In other words, the next order of
neural interactions in grConv must wait for the previous order’s to complete.
Therefore, a series of convolution operators are forced to be serialized,
which largely compromises the computing parallelism performance and
exacerbates computing latency, especially with GPUs.

Although such excessive intermediate states and data dependency bring new
challenges to high-order interactive neural networks, the above motivation
analysis also demystifies the core cognition mechanism of Transformer-like
models and the computation problems. To achieve a revolutionary design with
both cognition and computation performance in this work, we will show that
quadratic neurons can play an important role.

Figure 2: XOR Representation Capacity of Different Neurons
(colors: decision values, crosses/dots: positive/negative samples)

Figure 3: A Bottom-Up Overview of QuadraNet (i.e., from neuron level to model
level to system implementation level)

###  II-C Quadratic Neuron

A quadratic neuron brings more opportunities to discover nonlinear relations
among data dimensions because it performs a high-order interaction within the
neuron for its quadratic computation. Therefore, with a shallower network, a
quadratic neural network can learn a smoother and richer nonlinear
representation than conventional neural networks. For example, for a
generalized XOR classification problem, as shown in Fig. 2, a single-layer
perceptron with a quadratic neuron can solve it alone and gives even more
precise classification boundaries than a two-layer perceptron with ten hidden
neurons.

Neuron Format: In contrast to the conventional one-order neuron represented by
an inner product of input vector $\mathbf{x}\in\mathbb{R}^{n\times 1}$ and
weight vector $W\in\mathbb{R}^{1\times n}$ plus a bias $b$ and activation, the
quadratic neuron can be written as a quadratic polynomial of input
$\mathbf{x}$ in the formulation:

$\mathbf{y}=\mathbf{x}^{T}W_{q}\mathbf{x}+W_{c}\mathbf{x}+b,$ (4)

where $W_{q}$ and $W_{c}$ denote the quadratic- and linear-term weights,
respectively. In Eq. (4), $\mathbf{x}^{T}W_{q}$ can be regarded as an input-
dependent weight of $\mathbf{x}$; thus, it is evident that high-order
interaction is directly embedded in a quadratic neuron [13, 12]. This refines
the granularity of high-order interactions and their corresponding computation
from the model architecture level to the neuron level, allowing more flexible
adaption of high-order interactions. Therefore, we can leverage the neuron-
level high-order interaction as described in Eq. (4) to build high-order
neural networks with the mechanism in Eq. (3).

Opportunities and Challenges: By embedding high-order interactions into the
neuron level, ordinary quadratic neurons also generate a quadratic parameter
space as large as ($W_{q}\in\mathbb{R}^{n\times n}$), and its parameter volume
rises to $O$($n^{2}$+$n$) compared with $O$($n$) with conventional neuron-
based model implementations, since we now get an extra quadratic term weight
matrix with $n^{2}$ learnable parameters. Subsequently, the computing workload
in terms of Multiply Accumulate (MAC) operations rises to $O$($n^{2}$+$2n$).
Therefore, an inappropriate implementation of quadratic neural networks would
still cause certain computing issues, such as a significant parameter space
dilation and thus an increase in the amount of computation.

This paper undertakes a meticulous investigation of the quadratic neural
network, spanning from individual neurons to structural blocks, and ultimately
to the overall model design and automatic search. The primary aim is to
comprehensively address and optimize the parameter space dilation concern,
while simultaneously preserving the quadratic neurons’ robust high-order
neural interaction capabilities.

###  II-D Hardware-Aware Neural Architecture Search

Neural Architecture Search (NAS) has been extensively used in neural network
design and optimization. Particularly, hardware-aware NAS considers computing
cost, memory consumption, runtime latency, or power consumption as one of the
optimization objectives in model architecture search algorithms, thus
synthesizing a model that is well-balanced between cognition performance and
computation efficiency.

Earlier on, little attention from the NAS area is given to high-order
interactive neural networks due to the complex search space and large
evaluation cost. Only very recently, the popularity of Transformer-like models
triggered some works to build a NAS synthesized Vision Transformers [14].
However, these works only focused on cognition performance without considering
hardware awareness, not to mention the breakthrough in high-order neuron
designs.

In this work, we target to incorporate hardware awareness into the proposed
QuadraNet design and expand its generalization capacity. Therefore, we will
adopt mature NAS techniques into our design methodology to search for an
adaptive high-order interactive neural network architecture for various
deployment scenarios.

##  III QuadraNet Design Methodology

###  III-A Design Methodology Overview

Fig. 3 illustrates the overview of our proposed methodology, which consists of
the following major tasks:

(a) A novel computing-friendly quadratic neuron design is proposed. With
factorized low-rank quadratic terms, the parameter volumes and corresponding
workload could be significantly optimized at the neuron level first.

(b) Based on quadratic neurons, a quadratic convolution mechanism is built for
high-order interactions. Different from conventional approaches (e.g.,
$g^{r}Conv$), this mechanism could avoid large amounts of intermediate state
spaces and recursive data dependency, further retaining the computation
efficiency.

(c) QuadraBlock is designed as the basic high-order neural network structural
block. QuadraBlock achieves high-order spatial neural interaction with
quadratic depth-wise convolution and effectively unifies quadratic neurons
within interactive block for optimal cognition capacity.

(d) A high-order neural network architecture QuadraNet comprehensively
integrates the proposed techniques above and leads a design revolution for
Transformer-like models with outstanding cognition and computation
performance.

(e) Hardware-aware NAS is further incorporated to improve the architecture
generalization capacity of QuadraNet, namely, the adaptability to various
hardware computing constraints.

Eventually, external system implementation support is also well considered in
this work, such as a full-stack compiler and library integration to target
practical applications.

###  III-B Efficient Quadratic Neuron and Convolution Design

Quadratic Neuron: As defined in Eq. (4), the quadratic term of a conventional
quadratic neuron has a parameter volume of $n^{2}$, which limits its usage in
terms of both algorithm complexity and implementation difficulty. To address
this issue, instead of using a full-rank quadratic weight matrix $W_{q}$, we
propose to build the one-rank weight vectors $W_{a}$ and
$W_{b}\in\mathbb{R}^{1\times n}$ as shown in Fig. 3 (a), which can be seen as
a low-rank factorization of the quadratic weight term. Thus, the new quadratic
neuron can be reformulated as:

$\mathbf{y}=\mathbf{x}^{T}W_{a}^{T}W_{b}\mathbf{x}+W_{c}\mathbf{x}+b\\\
=(W_{a}\mathbf{x})^{T}(W_{b}\mathbf{x})+W_{c}\mathbf{x}+b.$ (5)

Such a factorization reduces the parameter space of quadratic neurons from
$O$($n^{2}$+$n$) to $O(3n)$ and the computation complexity from
$O$($n^{2}$+2$n$) to $O$(4$n$).

High-Order Interactive Quadratic Convolution: Because quadratic neurons have
powerful cognitive capacity and generalization ability, we can directly bring
them into a simple way of information extraction. In this work, we directly
project optimized quadratic neurons to the widely-used convolution operation
to achieve efficient high-order interactions. Specifically, we map the three
groups of weights $W_{a}$, $W_{b}$, and $W_{c}$ from a single quadratic neuron
to three convolutional filters and conduct a polynomial computation like in
Eq. (5). With such a mapping, we obtain a quadratic convolution operator.
Considering a mostly used depthwise separable convolution, the quadratic
computation of the $i^{th}$ input data of the $c^{th}$ channel is formulated
as follows:

$\displaystyle y_{Quad}^{(i,c)}$
$\displaystyle=\sum_{j,k\in\Omega_{i}}\sum_{c^{\prime}=1}^{C}(\omega_{a,i\rightarrow
j}T_{a}^{(c,c^{\prime})}x^{(j,c^{\prime})})(\omega_{b,i\rightarrow
k}T_{b}^{(c,c^{\prime})}x^{(k,c^{\prime})})$ (6)
$\displaystyle=\sum_{j\in\Omega_{i}}\sum_{c^{\prime}=1}^{C}\omega_{a,i\rightarrow
j}\sum_{k\in\Omega_{i}}\omega_{b,i\rightarrow
k}x^{(k,c^{\prime})}T_{q}^{(c,c^{\prime})}x^{(j,c^{\prime})}$
$\displaystyle=\sum_{j\in\Omega_{i}}\sum_{c^{\prime}=1}^{C}q(x)_{ij}T_{q}^{(c,c^{\prime})}x^{(j,c^{\prime})},$

where $\Omega_{i}$ denotes the convolution receptive field in the $i^{th}$
position, and $q(x)_{ij}=\omega_{a,i\rightarrow
j}\sum_{k\in\Omega_{i}}\omega_{b,i\rightarrow k}x^{(k,c^{\prime})}$ implicates
an input-adaptive weight. When multiplying $x$, it captures the high-order
spatial interaction. Notice that Eq. (6) forms a unified expression for high-
order interaction that well matches with Eq. (3). This indicates that the
quadratic convolution can capture high-order interactions, and therefore can
be seen as an approximation computation of self-attention and grConv.

Efficient High-Order Computation: Through the neuron- and convolution-level
design, the primary high-order quadratic computation mechanism could be
determined microscopically. From a tensor computation perspective, the
computing of this mechanism could be formulated as follows:

$X_{out}=f_{a}(X_{in})\odot f_{b}(X_{in})+f_{c}(X_{in}),$ (7)

where $\odot$ denotes the Hadamard product; $f_{a},f_{b},f_{c}$ denote
convolution operators (with weights of $W_{a}$, $W_{b}$, $W_{c}$); $X_{in}$
and $X_{out}$ denote the input and output tensors, respectively. Different
from [9], such a mechanism goes through a straightforward computation flow
without a recursive logic, and thus eliminating the aforementioned data
dependency and corresponding computing and memory costs.

As for the problem of intermediate states, we need to examine the ones
generated in Eq. (7). During forward propagation, four intermediate states are
generated — $f_{a}(X_{in})$,$f_{b}(X_{in})$,$f_{a}(X_{in})\odot
f_{b}(X_{in})$, and $f_{c}(X_{in})$, all of which require the same data volume
and 4$HWC$ intermediate state workload in total. The corresponding workload of
intermediate states in a quadratic convolution will therefore be $4\times$
compared with a conventional convolution operator (i.e., $HWC$) for better
cognition performance, but be much more efficient compared with the self-
attention mechanism (i.e., ($HW$+3)$\times$) for less computing costs.

Meanwhile, during backpropagation, the $4\times$ intermediate states can be
further reduced to $2\times$, since two intermediate states of
$f_{a}(X_{in})\odot f_{b}(X_{in})$ and $f_{c}(X_{in})$ can be immediately
released during computing before causing dominant workload. For the rest two
states — $f_{a}(X_{in})$ and $f_{b}(X_{in})$, weight updates of $W_{a}$ (or
$W_{b}$) only rely on the intermediate state of $f_{b}(X_{in})$ (or
$f_{a}(X_{in})$), because the gradient of the weight $W_{a}^{l}$ (or
$W_{b}^{l}$) can be computed by:

$\displaystyle\frac{\partial\mathscr{L}}{\partial W_{a}^{l}}=$
$\displaystyle\frac{\partial\mathscr{L}}{\partial X^{l+1}}\cdot\frac{\partial
X^{l+1}}{\partial(W_{a}^{l}X^{l})(W_{b}^{l}X^{l})}\cdot\frac{\partial(W_{a}^{l}X^{l})(W_{b}^{l}X^{l})}{\partial(W_{a}^{l}X^{l})}\cdot\frac{\partial(W_{a}^{l}X^{l})}{\partial
W_{a}^{l}}$ (8) $\displaystyle=$
$\displaystyle\frac{\partial\mathscr{L}}{\partial
X^{l+1}}\cdot(W_{b}^{l}X^{l})\cdot X^{l},$

where $X^{l}$ denotes the $l^{th}$ layer output. Moreover, the weight update
of the linear term $W_{c}$ does not rely on any intermediate states.
Therefore, only two intermediate states — $f_{a}(X_{in})$ and $f_{b}(X_{in})$
are actually counted for memory and computing costs.

In a word, quadratic convolution enables more efficient higher-order
computation by avoiding a recursive design for fewer data dependencies and
further optimizations for intermediate states. Quadratic convolution also
offers the possibility of external acceleration by simply mapping quadratic
neurons to the convolution operators, allowing us to directly leverage
existing compilers [15] and the quadratic computation library [12] to optimize
the computation.

###  III-C Efficient and Powerful QuadraBlock

QuadraBlock Design: Our QuadraBlock design is as shown in Fig. 3 (c). It is
divided into a spatial feature mixing layer and an inter-channel information
mixing layer: the spatial feature mixing layer is composed of a layer
normalization operator and a quadratic depth-wise convolution; the inter-
channel information mixing layer comprises a layer normalization operator, a
1$\times$1 convolution operator with a channel expansion coefficient $R$, a
GELU activation, and another 1$\times$1 convolution operator restoring the
number of channels.

Such a design achieves high computing efficiency by using only quadratic
convolutions in spatial information mixing. Because in the modern CNNs design,
the spatial feature mixing layers take less than 5% of the parameters and MACs
while performing important texture, outline, and shape depicting [16].
Therefore, inserting the quadratic convolutions into the spatial information
mixture has little effect on the overall number of parameters and computation
cost.

QuadraBlock for Cognition Capacity: QuadraBlock also unifies multiple
quadratic neurons using channle-wise information transportation and
activation. It enhances the cognition ability of quadratic neurons by
approximating a high-rank quadratic weight matrix within the QuadraBlock. This
enhancement approximation is achieved by leveraging the nonlinear activation
function to construct a linearly independent quadratic term weights matrix as
in Eq. (9). The quadratic term of the combination of quadratic neuron layers,
1$\times$1 convolution with expansion coefficient $R$ and activation function
$\sigma(\cdot)$ can be formulated as follows:

$\displaystyle y_{comb}^{c}$
$\displaystyle=\sum_{c^{\prime}=1}^{RC}\sigma(T_{q}^{(c,c^{\prime})}\mathbf{x}^{T}W_{a}^{T}W_{b}\mathbf{x})$
(9)
$\displaystyle\approx\sum_{c^{\prime}=1}^{C}T_{q}^{(c,c^{\prime})}\mathbf{x}^{T}(\sum_{r=1}^{R}W_{ar}^{T}W_{br})\mathbf{x},$

where
$(\sum_{r=1}^{R}W_{ar}^{T}W_{br})=W_{A}^{T}W_{B},W_{A},W_{B}\in\mathbb{R}^{n*R}$.

Thus, the delicately-designed QuadraBlock enhances the representation capacity
of quadratic weight matrix $W_{q}$ with $rank$=1 to a similar level as a
$rank$=$R$, retaining both low-rank computation efficiency and cognition
capability.

###  III-D QuadraNet Architecture Design:

Figure 3 (d) illustrates the macro design of QuadraNet. We follow full-fledged
four-stage pyramid architecture as [4, 16]. which contains 4-stage with
increasing channels and decreasing resolutions. In each stage, STEM [16] is
first utilized for downsampling while increasing the number of channels, then
the QuadraBlock is repeated for $n_{i}$ times in stage i. And following
configurations will be considered:

TABLE I: Cognition and Computation Performance Evaluation Model | #Para. | FLOPs | #Layer | Top-1 Acc. | Throughput
---|---|---|---|---|---
| (M) | (G) |  | (%) | (img/s)
Swin-T[4] | 29 | 4.5 | 12 | 81.3 | 1321
DeiT-S[17] | 22 | 4.6 | 12 | 79.8 | 1621
TNT-S[18] | 23.8 | 5.2 | 12 | 81.5 | 657
PVT-S[19] | 24.5 | 3.8 | 16 | 79.8 | 1433
T2T-ViT-14[20] | 21.5 | 4.8 | 14 | 81.5 | 1376
ConvNeXt-T[16] | 28 | 4.5 | 18 | 82.1 | 1944
HorNet-T[9] | 22 | 4 | 25 | 82.8 | 1254
QuadraNet25-T | 16.2 | 2.9 | 25 | 81.2 | 2433
QuadraNet36-T | 23.6 | 4.1 | 36 | 82.2 | 1971
Swin-S[4] | 50 | 8.7 | 24 | 83.0 | 827
ConvNeXt-S[16] | 50 | 8.7 | 36 | 83.1 | 1275
HorNet-S[9] | 50 | 8.8 | 25 | 83.8 | 813
QuadraNet36-S | 50.2 | 8.9 | 36 | 83.8 | 1117
Swin-B[4] | 88 | 15.4 | 24 | 83.5 | 662
ConvNeXt-B[16] | 89 | 15.4 | 36 | 83.8 | 969
HorNet-B[9] | 87 | 15.6 | 25 | 84.2 | 641
QuadraNet25-B | 61.2 | 11.1 | 25 | 83.6 | 1138
QuadraNet36-B | 90.4 | 15.8 | 36 | 84.1 | 892

Receptive field: We set the receptive field regarding convolution kernel size
to 7$\times$7 for a quadratic depth-wise convolution. Compared with
conventional 3$\times$3 practice, a large receptive field helps to capture
high-order spatial interactions.

Number of Channels: In each stage, we set the number of channels in a
gradually expanding way to [C, 2C, 4C, 8C] in each stage following common
practice, specifically, C=[16, 32, 64, 96, 128] to construct QuadraNet
L-XXS/XS/T/S/B models (for different layers and scales).

Number of Layers: To fairly compare it with baselines [4, 16, 9], we can set
the number of blocks as [3,3,27,3] and [2,3,18,2] in each stage to construct
QuadraNet 36/25-C with similar model volumes and number of layers.

TABLE II: Memory Consumption Evaluation Block | High-order | Traning Memory(MB)
---|---|---
|  | C=64 | C=128
SkipBlock(No DW-Conv) | $\times$ | 5515 | 12027
ConvBlock(Naive DW-Conv) | $\times$ | 7364 | 15711
SwinBlock[4] | ✓ | 12674 | 23472
HorBlock[9] | ✓ | 10724 | 20715
QuadraBlock | ✓ | 8108 | 17191
QuadraBlock*(w/ QuadraLib) | ✓ | 7794 | 16730

###  III-E Hardware-Aware Quadratic Neural Architecture Search

We further incorporate state-of-the-art hardware-aware NAS techniques to
further enhance QuadraNet’s architecture generalization capability and
hardware adaptability.

Search Space: Following a common NAS development approach, QuadraNet first
requires a new parameter search space. When we generalize the above
configuration options into automation search algorithm, we can therefore
obtain a unique high-order interactive search space. Following the original
ProxylessNAS [21], we specify a candidate set of QuadraBlock with different
kernel sizes {3, 5, 7} and expansion coefficients {2, 4}. We also add an
identity block that does nothing but provides a residual (skip) connection to
the candidate set to search for a shallow network.

Hardware Awareness: As demonstrated in previous sections, we have analyzed the
QuadraNet’s memory and computing costs, as well as optimization factors. Given
different systems’ hardware profiling regarding different neural network
operations, we can project different QuadraNet model designs into particular
resource consumption estimations. Therefore, the hardware-aware search
feedback can be well established in the NAS process considering the hardware
constraints.

Search Method: Leveraging these efforts, we modify ProxylessNAS [21] into a
QuadraNet-specified NAS framework by adapting the Search Space and feedback
configuration. Firstly, we train a supernet in defined QuadraNet Search Space.
Then with a gradient-based search strategy, considering hardware latency as a
part of the loss function. And eventually, we train the branch weights to
search for an appropriate architecture for the target device with particular
hardware constraints.

##  IV Experiment

Since current high-order interactive neural networks are primarily tailored
for image classification, we also apply the QuadraNet methodology to this
specific task. Specifically, we conduct experiments on ImageNet-1K dataset
[22], which contains 1.28M training images and 50K validation images from 1000
classes. We report ImageNet-1K’s Top-1 accuracy on the validation set. In the
following, we compare the proposed QuadraNet with state-of-the-art models
regarding both cognition and computation performance.

###  IV-A QuadraNet Evaluation with Manually Designed Models

Model Initialization: We train QuadraNets for 300 epochs (20 of which are
warmup) with AdamW [23] optimizer and a learning rate of 4e-3. We use a weight
decay of 0.05 and the same data augmentations as in [9]. The resolution of the
input images is 224$\times$224\. To stabilize the training, we apply gradient
clipping that limits the maximum value of gradient to 5, instead of the
superabundant batch normalization as in [12], which will hurt the
generalization ability of networks. We use the “Channel Last” memory layout
following [16].

Cognition Performance and Computing Costs: Table I evaluates the performance
of several manually designed QuadraNet models with other state-of-the-art
models at comparable FLOPs and parameters. Our models demonstrate competitive
performance across both Transformer-like and CNN-like architectures. With such
a cognition performance baseline, we measure their inference throughputs on an
A100 GPU. And results show that, with a similar image classification accuracy,
QuadraNets achieve 1.35$\times$ to 1.57$\times$ higher inference throughputs
than high-order works [4, 9]. Compared with ConvNeXt [16], QuadraNets achieve
higher cognitive performance and even higher throughputs most of the time.

Memory Consumption Analysis: We further compare another set of models with a
similar configuration except for different types of computation blocks to
evaluate QuadraNet’s memory consumption with 64 batch size, 24 layers, and
64/128 channels. Two low performant non-high-order models (SkipBlock and
ConvBlock) are also involved in examining the trade-off between high-order
interactive performance and costs. To make a fair comparison, the SkipBlock
completely removes the spatial neural interaction and related residual
connection, and the ConvBlock uses the same block design but only with the
conventional depth-wise convolution.

Table II shows that, under the same configuration, QuadraBlock achieves high-
order information interaction with a relatively smaller extra memory footprint
than SwinBlock and HorBlock. When coupled with accelerate library QuadraLib
[12] that further optimizes the intermediate state, QuadraBlock results in an
even smaller memory footprint.

TABLE III: Hardware-Aware QuadraNet NAS Searched Model | FLOPs | CPU | VPU | GPU | Top-1 Acc.
---|---|---|---|---|---
Proxyless-CPU | 333M | 299.6ms | 53.1ms | 1.8ms | 72.4
Proxyless-VPU | 275M | 379.1ms | 29.9ms | 1.7ms | 72.7
Proxyless-GPU | 1065M | 1521.1ms | 146.1ms | 2.9ms | 74.4
QuadraNet_CPU | 279M | 299.8ms | 34.7ms | 2.1ms | 73.1
QuadraNet_VPU | 312M | 437.1ms | 29.9ms | 1.9ms | 73.2
QuadraNet_GPU | 948M | 1334.2ms | 117.5ms | 2.9ms | 76.7

###  IV-B Hardware-Aware QuadraNet Model Generalization

We conduct hardware-aware NAS with the help of nn-Meter [24] for its powerful
latency estimation for different devices. We set the maximum latency for
mobile CPU (CortexA76), VPU (Intel Myriad VPU), and GPU (A100) to 300ms, 30ms,
and 3ms, respectively. Therefore, we compare the NAS models from regular
ProxylessNAS [21] and hardware-aware QuadraNet. The results are shown in Table
III. We observe that, given different hardware constraints, QuadraNet models
always perform better in terms of accuracy, FLOPs, and computation time. Fig 4
also gives an example of a searched model with less than 10ms latency on VPU.

Figure 4: Visualization of Searched QuadraNet_VPU_10ms Figure 5: Outstanding
Quadratic Performance w/ Various Configurations

###  IV-C Ablation Study

Superiority of Quadratic Neuron: We further quantitatively study the cognitive
capacity of quadratic neurons. Fig. 5 shows that, with different numbers of
channels and layers, the additional parameters introduced by quadratic neurons
can be well reflected in the consistent cognitive performance improvement,
which further demonstrates the effectiveness of higher-order spatial neural
interactions.

Receptive Field’s Impact: As high-order interaction is the essential mechanism
of this work, we would take the receptive field as a representative factor to
examine its correlation with QuadraNet design. We train the QuadraNet36-T with
different kernel sizes on 4 A100 GPUs to investigate the impact of the
receptive field. As shown in Table IV, larger receptive fields indeed improve
accuracy, but at the cost of considerable training time due to the inflated
computational scale and unoptimized operators. Therefore, this observation
provides a new perspective on the cognition and computation trade-off. As
mentioned above, we utilized a balance 7$\times$7 kernel, while more
investigation would be performed in the future.

Channel/Spatial-Wise Quadratic Neuron Utilization: After all, quadratic
neurons have more computation costs than ordinary ones. Thus appropriate
utilization of them becomes critical. In conventional model designs, spatial-
wise neurons only account for 5% while channel-wise accounts for 95%. When it
comes to quadratic settings, a dedicated analysis is shown in Table V to
identify which approach is more quadratic neuron-friendly. We can see that,
although spatial-wise utilization is very small, it has even better
effectiveness overall, highlighting a  QuadraNet design principle.

##  V Conclusion

In this work, we resented QuadraNet — a revolutionary neural network design
methodology that incorporates efficient and powerful high-order neural
interaction interactions. QuadraNet avoids the existing high-order
computational issues of intermediate states and data dependency down to the
neuron design level. Based on the design of efficient quadratic neurons into
QuadraBlock, we construct a new family of generic backbone networks
QuadraNets, and provide a new search space for high-order interactions for
hardware-aware NAS. Experiments demonstrate the superior cognitive capacity,
high throughputs, and low training memory footprints of the proposed
QuadraNet. Hardware-aware NAS further identifies better quadratic neural
network architecture under this new search space with quadratic neurons.

TABLE IV: Receptive Field’s Impact Receptive Field | 3*3 | 5*5 | 7*7 | 9*9 | 21*21 |  Global Filter [25]
---|---|---|---|---|---|---
Top-1 Acc. | 80.6 | 81.3 | 82.2 | 82.2 | 82.8 | 82.6
Training hours | 81.1 | 86.9 | 94.4 | 123.3 | 390.0 | 336.1
TABLE V: Channel-wise Quadratic Operator’s Impact  Models | #Params | #FLOPs | Top-1 Acc.
---|---|---|---
QuadraNet36-T | 23.6M | 4.1G | 82.2
+Quadratic 1*1 Conv1 | 44.6M | 7.8G | 82.3
+Quadratic 1*1 Conv2 | 65.6M | 11.5G | 82.7
QuadraNet36-B | 90.4M | 15.8G | 84.1
+Quadratic 1*1 Conv1 | 174.3M | 30.6G | 83.9
+Quadratic 1*1 Conv2 | 258.2M | 45.4G | 84.2

## References

  * [1] A. Krizhevsky, _et al._ , “ImageNet Classification with Deep Convolutional Neural Networks,” in _NeurIPS_ , 2012.
  * [2] K. He, _et al._ , “Deep Residual Learning for Image Recognition,” in _CVPR_ , 2016.
  * [3] A. Dosovitskiy, _et al._ , “An Image Is Worth 16x16 Words: Transformers for Image Recognition at Scale,” in _ICLR_ , 2021.
  * [4] Z. Liu, _et al._ , “Swin Transformer: Hierarchical Vision Transformer Using Shifted Windows,” in _ICCV_ , 2021.
  * [5] S. Mehta _et al._ , “MobileViT: Light-Weight, General-Purpose, and Mobile-Friendly Vision Transformer,” _arXiv:2110.02178_ , 2021.
  * [6] M. Sandler, _et al._ , “Mobilenetv2: Inverted Residuals and Linear Bottlenecks,” in _CVPR_ , 2018.
  * [7] F. D. Keles, _et al._ , “On the Computational Complexity of Self-Attention,” _arXiv:2209.04881_ , 2022.
  * [8] S. Wang, _et al._ , “Linformer: Self-attention with linear complexity,” _arXiv preprint arXiv:2006.04768_ , 2020.
  * [9] Y. Rao, _et al._ , “HorNet: Efficient High-Order Spatial Interactions with Recursive Gated Convolutions,” _arXiv:2207.14284_ , 2022.
  * [10] G. G. Chrysos, _et al._ , “P-Nets: Deep Polynomial Neural Networks,” in _CVPR_ , 2020.
  * [11] F. Fan, _et al._ , “A New Type of Neurons for Machine Learning,” _International Journal for Numerical Methods in Biomedical Engineering_ , 2018\.
  * [12] Z. Xu, _et al._ , “QuadraLib: A Performant Quadratic Neural Network Library for Architecture Optimization and Design Exploration,” in _MLSys_ , 2022.
  * [13] P. Mantini _et al._ , “CQNN: Convolutional Quadratic Neural Networks,” in _ICPR_ , 2021.
  * [14] C. Gong, _et al._ , “NASViT: Neural Architecture Search for Efficient Vision Transformers with Gradient Conflict-Aware SuperNet Training,” in _ICLR_ , 2021.
  * [15] T. Chen, _et al._ , “TVM: An Automated End-to-End Optimizing Compiler for Deep Learning,” in _OSDI_ , 2018.
  * [16] Z. Liu, _et al._ , “A ConvNet for the 2020s,” in _CVPR_ , 2022.
  * [17] H. Touvron, _et al._ , “Training Data-Efficient Image Transformers & Distillation through Attention,” in _ICML_ , 2021.
  * [18] K. Han, _et al._ , “Transformer in Transformer,” in _NeurIPS_ , 2021\.
  * [19] W. Wang, _et al._ , “Pyramid vision transformer: A versatile backbone for dense prediction without convolutions,” in _ICCV_ , 2021.
  * [20] L. Yuan, _et al._ , “Tokens-to-Token ViT: Training Vision Transformers from Scratch on ImageNet,” in _ICCV_ , 2021.
  * [21] H. Cai, _et al._ , “Proxylessnas: Direct neural architecture search on target task and hardware,” in _ICLR_ , 2018.
  * [22] J. Deng, _et al._ , “Imagenet: A large-scale hierarchical image database,” in _CVPR_ , 2009.
  * [23] I. Loshchilov _et al._ , “Decoupled weight decay regularization,” in _ICLR_ , 2018.
  * [24] L. L. Zhang, _et al._ , “nn-meter: towards accurate latency prediction of deep-learning model inference on diverse edge devices,” in _MobiSys_ , 2021\.
  * [25] Y. Rao, _et al._ , “Global filter networks for image classification,” _NeurIPS_ , 2021.


* Belczynski et al. (2010) Belczynski K., Bulik T., Fryer C. L., Ruiter A., Valsecchi F., Vink J. S., Hurley J. R., 2010, ApJ, 714, 1217
  * Blaauw (1961) Blaauw A., 1961, Bull. Astron. Inst. Netherlands, 15, 265
  * Blaes et al. (2002) Blaes O., Lee M. H., Socrates A., 2002, ApJ, 578, 775
  * Bond et al. (1982) Bond J. R., Arnett W. D., Carr B. J., 1982, in Rees M. J., Stoneham R. J., eds, NATO Advanced Study Institute (ASI) Series C Vol. 90, Supernovae: A Survey of Current Research. pp 303–311
  * Bondi & Hoyle (1944) Bondi H., Hoyle F., 1944, MNRAS, 104, 273
  * Braudo et al. (2022) Braudo J., Bear E., Soker N., 2022, MNRAS, 510, 4242
  * Brookshaw & Tavani (1993) Brookshaw L., Tavani M., 1993, ApJ, 410, 719
  * Cantiello et al. (2007) Cantiello M., Yoon S. C., Langer N., Livio M., 2007, A&A, 465, L29
  * Claeys et al. (2014) Claeys J. S. W., Pols O. R., Izzard R. G., Vink J., Verbunt F. W. M., 2014, A&A, 563, A83
  * Comerford & Izzard (2020) Comerford T. A. F., Izzard R. G., 2020, MNRAS, 498, 2957
  * Dewi & Tauris (2000) Dewi J. D. M., Tauris T. M., 2000, A&A, 360, 1043
  * Dominik et al. (2012) Dominik M., Belczynski K., Fryer C., Holz D. E., Berti E., Bulik T., Mandel I., O’Shaughnessy R., 2012, ApJ, 759, 52
  * Dong et al. (2014) Dong S., Katz B., Socrates A., 2014, ApJ, 781, L5
  * Dorozsmai & Toonen (2022) Dorozsmai A., Toonen S., 2022, arXiv e-prints, p. arXiv:2207.08837
  * Duffell et al. (2020) Duffell P. C., D’Orazio D., Derdzinski A., Haiman Z., MacFadyen A., Rosen A. L., Zrake J., 2020, ApJ, 901, 25
  * Eggleton (1983) Eggleton P. P., 1983, ApJ, 268, 368
  * Eggleton & Kiseleva-Eggleton (2001) Eggleton P. P., Kiseleva-Eggleton L., 2001, ApJ, 562, 1012
  * Evans (2011) Evans N. R., 2011, Bulletin de la Societe Royale des Sciences de Liege, 80, 663
  * Fabrycky & Tremaine (2007) Fabrycky D., Tremaine S., 2007, ApJ, 669, 1298
  * Ford et al. (2004) Ford E. B., Kozinsky B., Rasio F. A., 2004, ApJ, 605, 966
  * Fowler & Hoyle (1964) Fowler W. A., Hoyle F., 1964, ApJS, 9, 201
  * Fragione & Loeb (2019) Fragione G., Loeb A., 2019, MNRAS, 486, 4443
  * Fraley (1968) Fraley G. S., 1968, Ap&SS, 2, 96
  * Fryer et al. (2012) Fryer C. L., Belczynski K., Wiktorowicz G., Dominik M., Kalogera V., Holz D. E., 2012, The Astrophysical Journal, 749, 91
  * Georgy et al. (2011) Georgy C., Meynet G., Maeder A., 2011, A&A, 527, A52
  * Ghodla et al. (2022) Ghodla S., Eldridge J. J., Stanway E. R., Stevance H. F., 2022, MNRAS, 518, 860
  * Glanz & Perets (2021a) Glanz H., Perets H. B., 2021a, MNRAS, 500, 1921
  * Glanz & Perets (2021b) Glanz H., Perets H. B., 2021b, MNRAS, 507, 2659
  * Hamann & Koesterke (1998) Hamann W. R., Koesterke L., 1998, A&A, 335, 1003
  * Hamann et al. (1995) Hamann W. R., Koesterke L., Wessolowski U., 1995, A&A, 299, 151
  * Hamers & Thompson (2019) Hamers A. S., Thompson T. A., 2019, ApJ, 883, 23
  * Hamers et al. (2013) Hamers A. S., Pols O. R., Claeys J. S. W., Nelemans G., 2013, MNRAS, 430, 2262
  * Hamers et al. (2021) Hamers A. S., Rantala A., Neunteufel P., Preece H., Vynatheya P., 2021, MNRAS, 502, 4479
  * Hamers et al. (2022a) Hamers A. S., Glanz H., Neunteufel P., 2022a, ApJS, 259, 25
  * Hamers et al. (2022b) Hamers A. S., Perets H. B., Thompson T. A., Neunteufel P., 2022b, ApJ, 925, 178
  * Harrington (1968) Harrington R. S., 1968, AJ, 73, 190
  * Heger & Woosley (2002) Heger A., Woosley S. E., 2002, ApJ, 567, 532
  * Heggie (1975) Heggie D. C., 1975, MNRAS, 173, 729
  * Hills (1983) Hills J. G., 1983, ApJ, 267, 322
  * Hjellming & Webbink (1987) Hjellming M. S., Webbink R. F., 1987, ApJ, 318, 794
  * Holman et al. (1997) Holman M., Touma J., Tremaine S., 1997, Nature, 386, 254
  * Huang (1956) Huang S. S., 1956, AJ, 61, 49
  * Huang (1963) Huang S.-S., 1963, ApJ, 138, 471
  * Hurley et al. (2000) Hurley J. R., Pols O. R., Tout C. A., 2000, MNRAS, 315, 543
  * Hurley et al. (2002) Hurley J. R., Tout C. A., Pols O. R., 2002, MNRAS, 329, 897
  * Iben & Tutukov (1999) Iben Icko J., Tutukov A. V., 1999, ApJ, 511, 324
  * Innanen et al. (1997) Innanen K. A., Zheng J. Q., Mikkola S., Valtonen M. J., 1997, AJ, 113, 1915
  * Ivanova et al. (2013) Ivanova N., et al., 2013, A&ARv, 21, 59
  * Kaaret et al. (2017) Kaaret P., Feng H., Roberts T. P., 2017, ARA&A, 55, 303
  * Katz & Dong (2012) Katz B., Dong S., 2012, arXiv e-prints, p. arXiv:1211.4584
  * King et al. (2000) King A. R., Taam R. E., Begelman M. C., 2000, ApJ, 530, L25
  * Kinoshita & Nakai (1999) Kinoshita H., Nakai H., 1999, Celestial Mechanics and Dynamical Astronomy, 75, 125
  * Kiseleva et al. (1994) Kiseleva L. G., Eggleton P. P., Orlov V. V., 1994, MNRAS, 270, 936
  * Klencki et al. (2020) Klencki J., Nelemans G., Istrate A. G., Pols O., 2020, Astronomy & Astrophysics, 638, A55
  * Kobulnicky et al. (2014) Kobulnicky H. A., et al., 2014, ApJS, 213, 34
  * Kozai (1962) Kozai Y., 1962, AJ, 67, 591
  * Kroupa (2001) Kroupa P., 2001, MNRAS, 322, 231
  * Kummer et al. (2023) Kummer F., Toonen S., de Koter A., 2023, The Main Evolutionary Pathways of Massive Hierarchical Triple Stars (arXiv:2306.09400)
  * Lai & Muñoz (2022) Lai D., Muñoz D. J., 2022, arXiv e-prints, p. arXiv:2211.00028
  * Langer (2012) Langer N., 2012, Annual Review of Astronomy and Astrophysics, 50, 107–164
  * Laplace et al. (2020) Laplace E., Götberg Y., de Mink S. E., Justham S., Farmer R., 2020, A&A, 637, A6
  * Leigh et al. (2020) Leigh N. W. C., Toonen S., Portegies Zwart S. F., Perna R., 2020, MNRAS, 496, 1819
  * Lidov (1962) Lidov M. L., 1962, Planet. Space Sci., 9, 719
  * Liu et al. (2015) Liu B., Muñoz D. J., Lai D., 2015, MNRAS, 447, 747
  * Lower et al. (2018) Lower M. E., Thrane E., Lasky P. D., Smith R., 2018, Phys. Rev. D, 98, 083028
  * Lubow & Shu (1975) Lubow S. H., Shu F. H., 1975, ApJ, 198, 383
  * Madau & Dickinson (2014) Madau P., Dickinson M., 2014, ARA&A, 52, 415
  * Madau & Fragos (2017) Madau P., Fragos T., 2017, The Astrophysical Journal, 840, 39
  * Maeder (1987) Maeder A., 1987, A&A, 178, 159
  * Mandel & de Mink (2016) Mandel I., de Mink S. E., 2016, MNRAS, 458, 2634
  * Marchant et al. (2016) Marchant P., Langer N., Podsiadlowski P., Tauris T. M., Moriya T. J., 2016, A&A, 588, A50
  * Marchant et al. (2019) Marchant P., Renzo M., Farmer R., Pappas K. M. W., Taam R. E., de Mink S. E., Kalogera V., 2019, ApJ, 882, 36
  * Marchant et al. (2021) Marchant P., Pappas K. M. W., Gallegos-Garcia M., Berry C. P. L., Taam R. E., Kalogera V., Podsiadlowski P., 2021, A&A, 650, A107
  * Mardling & Aarseth (2001) Mardling R. A., Aarseth S. J., 2001, MNRAS, 321, 398
  * Martinez et al. (2022) Martinez M. A. S., Rodriguez C. L., Fragione G., 2022, ApJ, 937, 78
  * Mazeh & Shaham (1979) Mazeh T., Shaham J., 1979, A&A, 77, 145
  * McKernan et al. (2019) McKernan B., et al., 2019, ApJ, 884, L50
  * McKernan et al. (2020) McKernan B., Ford K. E. S., O’Shaugnessy R., Wysocki D., 2020, MNRAS, 494, 1203
  * Michaely & Perets (2014) Michaely E., Perets H. B., 2014, ApJ, 794, 122
  * Miller & Hamilton (2002) Miller M. C., Hamilton D. P., 2002, ApJ, 576, 894
  * Misner et al. (1973) Misner C. W., Thorne K. S., Wheeler J. A., 1973, Gravitation
  * Moe & Di Stefano (2017) Moe M., Di Stefano R., 2017, ApJS, 230, 15
  * Moody et al. (2019) Moody M. S. L., Shi J.-M., Stone J. M., 2019, ApJ, 875, 66
  * Moreno Méndez et al. (2022) Moreno Méndez E., De Colle F., López Cámara D., Vigna-Gómez A., 2022, arXiv e-prints, p. arXiv:2207.03514
  * Muñoz et al. (2019) Muñoz D. J., Miranda R., Lai D., 2019, ApJ, 871, 84
  * Naoz (2016) Naoz S., 2016, ARA&A, 54, 441
  * Naoz et al. (2013) Naoz S., Farr W. M., Lithwick Y., Rasio F. A., Teyssandier J., 2013, MNRAS, 431, 2155
  * Nieuwenhuijzen & de Jager (1990) Nieuwenhuijzen H., de Jager C., 1990, A&A, 231, 134
  * Öpik (1924) Öpik E., 1924, Publications of the Tartu Astrofizica Observatory, 25, 1
  * Paczynski (1976) Paczynski B., 1976, Symposium - International Astronomical Union, 73, 75–80
  * Paxton et al. (2011) Paxton B., Bildsten L., Dotter A., Herwig F., Lesaffre P., Timmes F., 2011, ApJS, 192, 3
  * Perets & Fabrycky (2009) Perets H. B., Fabrycky D. C., 2009, ApJ, 697, 1048
  * Perets & Kratter (2012) Perets H. B., Kratter K. M., 2012, ApJ, 760, 99
  * Peters (1964) Peters P. C., 1964, Phys. Rev., 136, B1224
  * Petrovich (2015) Petrovich C., 2015, ApJ, 799, 27
  * Pijloo et al. (2012) Pijloo J. T., Caputo D. P., Portegies Zwart S. F., 2012, MNRAS, 424, 2914
  * Pols et al. (1998) Pols O. R., Schröder K.-P., Hurley J. R., Tout C. A., Eggleton P. P., 1998, MNRAS, 298, 525
  * Portegies Zwart & Verbunt (1996) Portegies Zwart S. F., Verbunt F., 1996, A&A, 309, 179
  * Rakavy & Shaviv (1967) Rakavy G., Shaviv G., 1967, ApJ, 148, 803
  * Reimers (1975) Reimers D., 1975, Memoires of the Societe Royale des Sciences de Liege, 8, 369
  * Renzo et al. (2020) Renzo M., Farmer R., Justham S., Götberg Y., de Mink S. E., Zapartas E., Marchant P., Smith N., 2020, A&A, 640, A56
  * Renzo et al. (2021) Renzo M., et al., 2021, ApJ, 919, 128
  * Riley et al. (2021) Riley J., Mandel I., Marchant P., Butler E., Nathaniel K., Neijssel C., Shortt S., Vigna-Gómez A., 2021, MNRAS, 505, 663
  * Rodriguez & Antonini (2018) Rodriguez C. L., Antonini F., 2018, ApJ, 863, 7
  * Rozner & Perets (2022) Rozner M., Perets H. B., 2022, ApJ, 931, 149
  * Sana et al. (2012) Sana H., et al., 2012, Science, 337, 444
  * Sana et al. (2014) Sana H., et al., 2014, ApJS, 215, 15
  * Sander & Vink (2020) Sander A. A. C., Vink J. S., 2020, MNRAS, 499, 873–892
  * Schrøder et al. (2021) Schrøder S. L., MacLeod M., Ramirez-Ruiz E., Mandel I., Fragos T., Loeb A., Everson R. W., 2021, arXiv e-prints, p. arXiv:2107.09675
  * Shappee & Thompson (2013) Shappee B. J., Thompson T. A., 2013, ApJ, 766, 64
  * Silsbee & Tremaine (2017) Silsbee K., Tremaine S., 2017, ApJ, 836, 39
  * Smeyers & Willems (2001) Smeyers P., Willems B., 2001, A&A, 373, 173
  * Soberman et al. (1997) Soberman G. E., Phinney E. S., van den Heuvel E. P. J., 1997, A&A, 327, 620
  * Soker & Bear (2021) Soker N., Bear E., 2021, MNRAS, 505, 4791
  * Song et al. (2016) Song H. F., Meynet G., Maeder A., Ekström S., Eggenberger P., 2016, A&A, 585, A120
  * Stegmann et al. (2022a) Stegmann J., Antonini F., Schneider F. R. N., Tiwari V., Chattopadhyay D., 2022a, Phys. Rev. D, 106, 023014
  * Stegmann et al. (2022b) Stegmann J., Antonini F., Moe M., 2022b, MNRAS, 516, 1406
  * Stevenson et al. (2019) Stevenson S., Sampson M., Powell J., Vigna-Gómez A., Neijssel C. J., Szécsi D., Mandel I., 2019, ApJ, 882, 121
  * Stone et al. (2017) Stone N. C., Metzger B. D., Haiman Z., 2017, MNRAS, 464, 946
  * Swaruba Rajamuthukumar et al. (2022) Swaruba Rajamuthukumar A., Hamers A., Neunteufel P., Pakmor R., de mink S. E., 2022, arXiv e-prints, p. arXiv:2211.04463
  * The LIGO Scientific Collaboration et al. (2021) The LIGO Scientific Collaboration et al., 2021, arXiv e-prints, p. arXiv:2111.03606
  * Thompson (2011) Thompson T. A., 2011, ApJ, 741, 82
  * Tiede et al. (2020) Tiede C., Zrake J., MacFadyen A., Haiman Z., 2020, ApJ, 900, 43
  * Tokovinin (2010) Tokovinin A., 2010, VizieR Online Data Catalog, p. J/MNRAS/389/925
  * Tokovinin (2014) Tokovinin A., 2014, AJ, 147, 87
  * Tokovinin et al. (2006) Tokovinin A., Thomas S., Sterzik M., Udry S., 2006, A&A, 450, 681
  * Toonen et al. (2012) Toonen S., Nelemans G., Portegies Zwart S., 2012, A&A, 546, A70
  * Toonen et al. (2016) Toonen S., Hamers A., Portegies Zwart S., 2016, Computational Astrophysics and Cosmology, 3, 6
  * Toonen et al. (2018) Toonen S., Perets H. B., Hamers A. S., 2018, A&A, 610, A22
  * Toonen et al. (2020) Toonen S., Portegies Zwart S., Hamers A. S., Bandopadhyay D., 2020, A&A, 640, A16
  * Toonen et al. (2022) Toonen S., Boekholt T. C. N., Portegies Zwart S., 2022, A&A, 661, A61
  * Tutukov & Yungelson (1979) Tutukov A., Yungelson L., 1979, in Conti P. S., De Loore C. W. H., eds, Vol. 83, Mass Loss and Evolution of O-Type Stars. pp 401–406
  * Ulrich & Burger (1976) Ulrich R. K., Burger H. L., 1976, ApJ, 206, 509
  * Vassiliadis & Wood (1993) Vassiliadis E., Wood P. R., 1993, ApJ, 413, 641
  * Verbunt et al. (2017) Verbunt F., Igoshev A., Cator E., 2017, Astronomy & Astrophysics, 608, A57
  * Vigna-Gómez et al. (2021) Vigna-Gómez A., Toonen S., Ramirez-Ruiz E., Leigh N. W. C., Riley J., Haster C.-J., 2021, ApJ, 907, L19
  * Vigna-Gómez et al. (2022) Vigna-Gómez A., Liu B., Aguilera-Dena D. R., Grishin E., Ramirez-Ruiz E., Soares-Furtado M., 2022, MNRAS, 515, L50
  * Vink & de Koter (2005) Vink J. S., de Koter A., 2005, A&A, 442, 587
  * Vink et al. (2001) Vink J. S., de Koter A., Lamers H. J. G. L. M., 2001, A&A, 369, 574
  * Wen (2003) Wen L., 2003, ApJ, 598, 419
  * Woosley (2017) Woosley S. E., 2017, ApJ, 836, 244
  * Woosley & Weaver (1982) Woosley S. E., Weaver T. A., 1982, in Rees M. J., Stoneham R. J., eds, NATO Advanced Study Institute (ASI) Series C Vol. 90, Supernovae: A Survey of Current Research. p. 79
  * Yoon & Langer (2005) Yoon S.-C., Langer N., 2005, A&A, 443, 643
  * Yoon et al. (2006) Yoon S.-C., Langer N., Norman C., 2006, A&A, 460, 199
  * Yoshida et al. (2016) Yoshida T., Umeda H., Maeda K., Ishii T., 2016, MNRAS, 457, 351
  * de Kool et al. (1987) de Kool M., van den Heuvel E. P. J., Pylyser E., 1987, A&A, 183, 47
  * de Mink & King (2017) de Mink S. E., King A., 2017, ApJ, 839, L7
  * de Mink & Mandel (2016) de Mink S. E., Mandel I., 2016, MNRAS, 460, 3545–3553
  * de Mink et al. (2009) de Mink S. E., Cantiello M., Langer N., Pols O. R., Brott I., Yoon S. C., 2009, A&A, 497, 243
  * de Vries et al. (2014) de Vries N., Portegies Zwart S., Figueira J., 2014, MNRAS, 438, 1909
  * du Buisson et al. (2020) du Buisson L., et al., 2020, MNRAS, 499, 5941
  * van den Heuvel et al. (2017) van den Heuvel E. P. J., Portegies Zwart S. F., de Mink S. E., 2017, MNRAS, 471, 4256
  * von Zeipel (1910) von Zeipel H., 1910, Astronomische Nachrichten, 183, 345

##  Appendix A Additional figures

In this section, we present the low metallicity model (i.e. Z = 0.0005)
counterparts of some of the figures presented in the main text. In Fig. 4, we
show the mass ratios at the onset of TMT and the amount of (relative) mass
transferred towards the inner binary. In Fig. 14, we show the cumulative
distribution of eccentricities for CHE triples that experience TMT at the
onset of the mass transfer phase. In 15, we show the $M_{\rm 2,ZAMS}$
distribution of TMT sources, distinguishing them based on the evolutionary
phase of the inner binary. In Fig. 16 we show the distribution of initial
inner pericenters of CHE triples, distinguishing systems based on the maxium
inner eccentricity reached during evolution (left panel), and the based on the
evolutionary channel (right panel). In Fig 17, we show the outer pericentre
before and after the TMT episode for systems with MS-MS inner binaries (upper
panel) and BH-BH inner binaries (lower panel) at the onset of the mass
transfer phase.

Figure 13: The same figure as Fig. 4 but for our model at Z = 0.0005.  Figure
14: The same figure as Fig. 5 but for our model at Z = 0.0005.  Figure 15: The
same as Fig. 7 but for our model at Z = 0.0005

Figure 16: The same as Fig. 8 but at Z = 0.0005 Figure 17: The orbital
separation of the outer orbit before and at the onset of TMT for MS-MS inner
(upper panel) and BH-BH inner binary accretors (lower panel) as calculated
with equation 30 at Z = 0.0005.

##  Appendix B Calculation of birth and event rates

Throughout the paper, we estimate the:

  1. 1.

Formation efficiency (equation 31)

  2. 2.

Birth rate density (equation 38)

  3. 3.

Merger rate density (equation 45)

for each identified evolutionary channels. In this section, we discuss in
detail how we determine these quantities.

(i) Formation efficiency: The formation efficiency expresses the number of
ZAMS stellar systems formed that will evolve according to a specific
evolutionary channel as a fraction of all ZAMS stellar systems formed. We
calculate this quantity as:

$\epsilon_{\rm formation}=f_{\rm pm}\cdot\frac{N_{\rm channel}}{N_{\rm
simulated}},$ (31)

where $N_{\rm channel}$ is the number simulated systems that evolves according
to the channel of interest, $N_{\rm\rm simulated}$ the total number of sampled
systems, and $f_{\rm pm}$ is the portion of the simulated parameter space with
respect to the complete parameter space, that is:

$f_{\rm pm}=f_{\rm triple}\cdot f_{\rm M_{\rm 1,ZAMS}}\cdot f_{\rm q,in}\cdot
f_{\rm q,out}\cdot f_{\rm a,in}\cdot f_{\rm a,out},$ (32)

where $f_{\rm triple}$ is the assumed triple fraction, $f_{\rm M_{\rm
1,ZAMS}}$ is the fraction of the simulated parameter space of primary masses:

$f_{\rm M_{\rm 1,ZAMS}}=\frac{\int_{20\,M_{\odot}}^{100\,M_{\odot}}M_{\rm
1,ZAMS}^{-2.3}dm}{\int_{0.08\,M_{\odot}}^{0.5\,M_{\odot}}M_{\rm
1,ZAMS}^{-1.3}dm+\int_{0.5\,M_{\odot}}^{100\,M_{\odot}}M_{\rm
1,ZAMS}^{-2.3}dm},$ (33)

where we assumed that the absolute minimum stellar mass is $M_{\rm
ZAMS,min}=0.08\,M_{\odot}$ and the absolute maximum stellar mass is $M_{\rm
ZAMS,max}=100\,M_{\odot}$, and as explained in section 2.6, we sample primary
masses in the range of 20-100$\,M_{\odot}$. The fraction of the simulated
parameter space of inner mass ratios is:

$f_{\rm q,in}=\frac{1.0-0.7}{1.0-0.0},$ (34)

since the distribution of (inner and outer) mass ratios is assumed to be
uniform. In equation 34, we assume that inner mass ratios of hierarchical
triples have an interval of (0,1] and we sample from the interval of [0.7,1].
The fraction of the simulated parameter space of outer mass ratios is

$f_{\rm q,out}=\frac{1.0-0.1}{1.0-0.0},$ (35)

where we assume that outer mass ratios triples have an interval of (0,1] and
we sample from the interval of [0.1,1]. The fraction of the simulated
parameter space of inner semimajor axis is:

$f_{\rm
a,in}=\frac{\rm{log}_{10}(40\,R_{\odot})-\rm{log}_{10}(14\,R_{\odot})}{\rm{log}_{10}(10^{5}\,R_{\odot})-\rm{log}_{10}(14\,R_{\odot})},$
(36)

since the distribution of (inner and outer) semimajor axis is assumed to be
uniform in a logarithmic space. We assume that inner mass semimajor axes of
all triples range from 14$\,R_{\odot}$ to $10^{5}\,R_{\odot}$ and we sample
from the interval of [14,40] $R_{\odot}$. Finally, the fraction of the
simulated parameter space of outer semimajor axis is:

$f_{\rm
a,out}=\frac{\rm{log}_{10}(10^{5}\,R_{\odot})-\rm{log}_{10}(10^{2}\,R_{\odot})}{\rm{log}_{10}(10^{5}\,R_{\odot})-\rm{log}_{10}(10^{2}\,R_{\odot})},$
(37)

where we assume that inner mass semimajor axes of all triples range from
$10^{2}\,R_{\odot}$ to $10^{5}\,R_{\odot}$ and we sample from the enitre
interval. Equation 32 for channels involving isolated binaries reduces to
$f_{\rm pm}=f_{\rm binary}\cdot f_{\rm M_{\rm 1,ZAMS}}\cdot f_{\rm q,in}\cdot
f_{\rm a,in}$.

ii) Birth rate density: The birth rate density gives the number density of
ZAMS stellar systems in the local universe (that is at redshift $z\approx 0$),
which will evolve according to a specific channel. We calculate the birth rate
of systems in a certain channel as:

$R_{\rm birth}=\sum_{Z_{i}}\frac{\textrm{SFRd}^{*}(Z_{i},z_{\rm
ZAMS}=0)}{\tilde{M}}\cdot\epsilon_{\rm formation},$ (38)

where we sum over the two metallicity values, at which we performed our
simulations; Z = 0.005 and Z = 0.0005. $\rm{SFRd}^{*}(Z,z)$ is defined as the
metallicity-specific star formation rate density, and it gives the stellar
mass formed within a metallicity range $Z_{\rm low}\leq Z\leq Z_{\rm high}$ at
redshift $z$:

$\rm{SFRd}^{*}(Z,z)=\int_{Z_{\rm low}}^{Z_{\rm high}}f_{\rm
met}(Z,z)\rm{SFRd}(z)dZ,$ (39)

where $Z_{\rm low}$ and $Z_{\rm high}$ are 0.0015 ($10^{-10}$) and 0.01
(0.0015), respectively, for our model with Z = 0.005 (Z = 0.0005). Here, Z =
0.0015 is the midpoint between Z = 0.005 and Z = 0.0005 in logarithmic space,
Z = 0.01 is the highest metallicity at which CHE binaries can still form GW
sources at appreciable numbers and Z = $10^{-10}$ is an arbitrarily chosen,
extremely low metallicity value. In equation 39, $\rm{SFRd}(z)$ is the star
formation rate density, and we use the model from Madau & Dickinson (2014):

$\textrm{SFRd}(z)=\frac{0.01\cdot(1+z)^{2.6}}{1+((1+z)/3.2)^{6.2}}\,M_{\odot}\rm{yr^{-1}}Mpc^{-3},$
(40)

and $f_{\rm met}(Z,z)$ is the metallicity distribution of the stellar mass
formed. This quantity is also redshift dependent and assumed to follow a log-
normal distribution (Madau & Fragos, 2017):

$f_{\rm
met}(Z,z)=\frac{1}{\sigma\sqrt{2\pi}}\exp\left(\frac{(\log_{10}(Z)-\mu(z))^{2}}{2\sigma^{2}}\right),$
(41)

with a standard deviation of $\sigma=0.5$ and with a redshift-dependent mean
metallicity $\mu(z)=\log_{10}(Z_{\odot}\cdot
10^{0.153-0.074z^{1.34}})-0.5\rm{ln}(10)\sigma^{2}$. Finally, the term
$\tilde{M}$ in equation 38 is the average mass of all stellar systems and we
calculate this as:

$\begin{split}&\tilde{M}=f_{\rm single}\cdot\tilde{M}_{\rm 1,ZAMS}+\\\ &f_{\rm
binary}\cdot\int_{0}^{1}(1+q_{\rm in})\tilde{M}_{\rm 1,ZAMS}dq_{\rm in}+\\\
&f_{\rm triple}\cdot\int_{0}^{1}\int_{0}^{1}(1+q_{\rm in})(1+q_{\rm
out})\tilde{M}_{\rm 1,ZAMS}dq_{\rm in}dq_{\rm out},\end{split}$ (42)

where we have defined $\tilde{M}_{\rm 1,ZAMS}$, as the average mass of the
primary, i.e.:

$\tilde{M}_{\rm 1,ZAMS}=\int_{0.08\,M_{\odot}}^{100\,M_{\odot}}M_{\rm
1,ZAMS}f_{\rm IMF}dM_{\rm 1,ZAMS}$ (43)

where $f_{\rm IMF}$ is the normalised, piecewise continuous initial mass
function of Kroupa (2001), $f_{\rm single}$ and $f_{\rm binary}$ are the
single and binary fractions, respectively. We neglect higher order systems,
such that $f_{\rm triple}=1-f_{\rm single}-f_{\rm binary}$.

We note that we also assume that the binary and triple fractions are
independent on the primary mass of the system (which is clearly not consistent
with observations, see e.g. Moe & Di Stefano, 2017, but commonly assumed in
population synthesis studies as a simplification). Assuming flat mass ratio
distributions for both the inner and outer binary, equation 42 becomes:

$\tilde{M}=\left(f_{\rm single}+\frac{3}{2}f_{\rm binary}+\frac{9}{4}\cdot
f_{\rm triple}\right)\cdot\tilde{M}_{\rm 1,ZAMS},$ (44)

The term $\textrm{SFRd}^{*}(Z_{i},z_{\rm ZAMS}=0)/\tilde{M}$ in equation 38
then gives the average number of stars formed at redshift $z=0$ in a
metallicity range of $Z_{\rm i,low}\leq Z\leq Z_{\rm i,high}$. Multiplying
this term with $\epsilon_{\rm formation}$ gives the number of systems formed
in a given formation channel as a fraction of all systems formed in the above
mentioned metallicity range for a given star formation history model. Summing
these values over all of our metallicity bins therefore yields the total birth
rate of systems in a specific channel.

iii) Merger rate density: The merger rate density gives the rate density of a
given astrophysical event (such as GW transients from coalescing double
compact objects) in the local universe. The main difference between the birth
and merger rate is due to the considerable delay time between the formation of
the stellar system and the occurrence of the GW merger. For example, if the
delay time for a GW source at $z=0$ is $t_{\rm delay}=$ 10.5 Gyrs, then the
redshift at ZAMS of its progenitor systems is $z_{\rm ZAMS}\approx 2$, at
which the star formation rate density is an order of magnitude higher with
respect to its value at $z=0$ (e.g. see models of Madau & Dickinson, 2014). We
determine the merger rate density at $z=0$ as:

$R_{\rm
event}=\sum_{Z_{i}}\int_{0\,\rm{Gyr}}^{13.5\,\rm{Gyr}}\frac{\textrm{SFRd}^{*}(Z_{i},z_{\rm
ZAMS}(t_{\rm delay}))}{\tilde{M}}\cdot\tilde{\epsilon}(t_{\rm delay})dt_{\rm
delay},$ (45)

where $z_{\rm ZAMS}$ is the redshift at which the progenitor of a given
astrophysical event is formed (and therefore it is a function of delay time),
$\tilde{\epsilon}$ is the number of astrophysical events occurring at z = 0
with a delay time of $t_{\rm delay}$ as a fraction of all ZAMS stellar systems
formed at z = $z_{\rm ZAMS}$. We determine $z_{\rm ZAMS}$ for a given delay
time via the standard relation for lookback time:

$t_{\rm delay}=\frac{1}{H_{0}}\int_{z=0}^{z_{\rm
ZAMS}}\frac{dz^{\prime}}{(1+z^{\prime})E(z^{\prime})}dz^{\prime},$ (46)

where $E(z)=\sqrt{\Omega_{m}(1+z)^{3}+\Omega_{\lambda}}$, with
$\Omega_{M}=0.3$, $\Omega_{\lambda}=0.7$ and
$H_{0}=70\rm{kms^{-1}}\rm{Mpc}^{-1}$.

We note that our merger rate density should be only considered as an order of
magnitude estimate at best. This imprecision is due to several uncertainties
in stellar physics and, notably, the limited density of our metallicity grid.
We performed simulations only at two metallicities to determine the merger
rate density. However, the formation efficiency and delay times of GW sources
originating from CHE systems is expected to be sensitively dependent on
metallicity.

In particular, we overestimate the delay times for GW sources formed at
$0.001<Z\leq 0.005$, which in turn leads to an overestimation of the merger
rate density at z = 0. This is because, we represent all systems formed in
this metallicity range with our models at Z = 0.005, at which the stellar
winds are stronger and therefore lead to wider BH-BH binaries. The longer time
delays imply that GW sources merging at z = 0 are predicted to have formed at
a larger redshift, at which the star formation rate is higher. In particular,
Madau & Dickinson 2014 predicts that the cosmic star formation rate
montonically increases up to $z\sim 2$. This could also explain why our merger
rate is a factor of two higher than predicted by (Riley et al., 2021).
Similarly, we underestimate the delay times for GW sources formed at
$0.0005<Z\leq 0.001$, and therefore we might underestimate the merger rate
densities for such systems. In particular, this could mean that the merger
rate density of the TMT with a MS-MS accretor channel (discussed in section
5.5) could be significantly higher than predicted (shown in Table 3).