# SubGen: Token Generation in Sublinear Time and Memory

Amir Zandieh111Equal Contributions
Independent Researcher
<EMAIL_ADDRESS>Insu Han∗
Yale University
<EMAIL_ADDRESS>Vahab Mirrokni
Google Research
<EMAIL_ADDRESS>Amin Karbasi
Yale University, Google Research
<EMAIL_ADDRESS>
###### Abstract

Despite the significant success of large language models (LLMs), their
extensive memory requirements pose challenges for deploying them in long-
context token generation. The substantial memory footprint of LLM decoders
arises from the necessity to store all previous tokens in the attention
module, a requirement imposed by key-value (KV) caching. In this work, our
focus is on developing an efficient compression technique for the KV cache.
Empirical evidence indicates a significant clustering tendency within key
embeddings in the attention module. Building on this key insight, we have
devised a novel caching method with sublinear complexity, employing online
clustering on key tokens and online $\ell_{2}$ sampling on values. The result
is a provably accurate and efficient attention decoding algorithm, termed
SubGen. Not only does this algorithm ensure a sublinear memory footprint and
sublinear time complexity, but we also establish a tight error bound for our
approach. Empirical evaluations on long-context question-answering tasks
demonstrate that SubGen significantly outperforms existing and state-of-the-
art KV cache compression methods in terms of performance and efficiency.

##  1 Introduction

Large Language Models (LLMs) [1, 23] play a crucial role in various natural
language processing applications, including dialog systems [22, 5], coding
assistance [4, 18], and image/video generations from text [17, 11]. All of
these models rely on the transformer architecture, with the attention
mechanism serving as the key component.

To fully harness the capabilities of LLMs, they must demonstrate both
efficiency and accuracy in generating long sequences. In practical
applications, deploying LLMs to generate tokens in an autoregressive manner
involves a sequential decoding process, where attention is dynamically applied
to each newly generated token. This process effectively constructs the output
sequence in a streaming manner, one token at a time. Therefore, as the
sequence grows, the model has to produce contextually relevant and coherent
content.

A common method for autoregressive attention decoding involves the use of key-
value (KV) caching, where key and value pairs from all preceding tokens are
cached and reused to prevent redundant computations. However, this approach
faces memory constraints, particularly when handling long sequences. In
particular, the memory requirements and runtime for generating each new token
increase linearly with context size, posing a significant challenge for
efficient processing of extensive sequences. This linear scaling directly
impedes practical applicability in real-world scenarios, such as chat systems,
where large contexts are often encountered.

In this work, we delve into the primary computational and memory bottleneck of
token generation. We propose SubGen, a novel approach designed to
significantly reduce the memory and runtime complexity of token generation,
moving from conventional linear growth to sublinear scale. To summarize, our
goal is to answer the following question:

_Can we approximate the attention output in decoding phase
in sublinear space/time complexity in context length?_

###  1.1 Related Work

Recent studies have underscored the need for efficient token generation,
particularly with the rise of long-range context datasets. Several recent
works have developed efficient strategies for compressing the KV cache. Zhang
et al. [28] proposed a greedy-type eviction algorithm that dynamically keeps
at most $k\ll n$ token embeddings based on the accumulated attention scores
where they refer to the Heavy Hitter Oracle (H2O). Liu et al. [14] empirically
observed that tokens with initially high attention scores tend to stay high
during the future generation process. Motivated by this observation, the
authors proposed a strategy that only keeps the most recent and pivotal tokens
whose attention scores are higher than a threshold. Ge et al. [9] proposed an
adaptive method of KV cache compression which identifies the intrinsic
structures of attention heads and uses them to determine the optimal
compression policy. Xiao et al. [26] observed that a simple eviction mechanism
that keeps only the first few and last few tokens does not degrade much the
decoding quality. They additionally proposed a fine-tuning method to solve
performance degradation from their method. Liu et al. [15] developed an
algorithm that reduces the generation latency by exploiting contextual
sparsity. In addition to algorithmic acceleration, there has also been a line
of work optimizing hardware resource configurations [20, 12]. However, to the
best of our knowledge, none of these works have achieved an efficient method
for KV cache with fully sublinear-time memory space.

On the lower bound side, achieving subquadratic amortized runtime for
producing output embeddings for $n$ tokens in the worst-case instances is
likely impossible without making assumptions about the input tokens [2, 19].
Therefore, to achieve fast runtime, it is necessary to rely on certain
assumptions about the input tokens.

###  1.2 Streaming Attention Problem

Deployment of LLMs involves performing attention decoding in a streaming
fashion. More precisely, the stream of tokens is represented as a sequence of
vector triplets
$({\bm{q}}_{1},{\bm{k}}_{1},{\bm{v}}_{1}),({\bm{q}}_{2},{\bm{k}}_{2},{\bm{v}}_{2}),\ldots({\bm{q}}_{n},{\bm{k}}_{n},{\bm{v}}_{n})$,
where ${\bm{q}}_{i},{\bm{k}}_{i},{\bm{v}}_{i}\in\mathbb{R}^{d}$ are queries,
keys, and values of the attention mechanism and $n$ is the total number of
tokens in the stream so far either in prompt or generation. 222We denote
vectors with lowercase boldface letters, e.g., ${\bm{v}}$, matrices with
uppercase boldface letters, e.g., ${\bm{M}}$, and sets with calligraphy
uppercase letters, e.g., $\mathcal{S}$. The operator norm of a matrix is
denoted as $\left\|\cdot\right\|_{{op}}$.  The objective of streaming
attention decoding is to compute the following:

$\mathrm{Attn}({\bm{q}}_{n},{\bm{K}}_{n},{\bm{V}}_{n})=\mathtt{softmax}({\bm{K}}_{n}\cdot{\bm{q}}_{n})^{\top}\cdot{\bm{V}}_{n},$
(1)

where ${\bm{K}}_{n},{\bm{V}}_{n}\in\mathbb{R}^{n\times d}$ are matrices
defined by stacking the keys and values in their respective rows:

${\bm{K}}_{n}:=\begin{bmatrix}{\bm{k}}_{1}^{\top}\\\ {\bm{k}}_{2}^{\top}\\\
\vdots\\\
{\bm{k}}_{n}^{\top}\end{bmatrix},~{}~{}~{}{\bm{V}}_{n}:=\begin{bmatrix}{\bm{v}}_{1}^{\top}\\\
{\bm{v}}_{2}^{\top}\\\ \vdots\\\ {\bm{v}}_{n}^{\top}\end{bmatrix}.$ (2)

The output $\mathrm{Attn}({\bm{q}}_{n},{\bm{K}}_{n},{\bm{V}}_{n})$ is then
used for predicting the next token and its token embedding is applied to a
transformer model and a new stream pair
$({\bm{q}}_{n+1},{\bm{k}}_{n+1},{\bm{v}}_{n+1})$ is generated. However,
storing these values and keys requires $O(nd)$ memory, posing a significant
space complexity challenge for long-context models with large $n$.

###  1.3 Overview of Contributions

We propose SubGen, an efficient method that accurately approximates the
attention decoder’s output in ?? while retaining only a small (sublinear)
subset of keys and values in the cache. In particular, SubGen computes an
estimator ${\bm{z}}_{n}$ for
$\mathrm{Attn}({\bm{q}}_{n},{\bm{K}}_{n},{\bm{V}}_{n})$ in sublinear time and
memory such that the error is bounded as follows:

$\displaystyle\left\|{\bm{z}}_{n}-\mathrm{Attn}({\bm{q}}_{n},{\bm{K}}_{n},{\bm{V}}_{n})\right\|_{2}\leq\varepsilon\left\|\mathtt{softmax}({\bm{K}}_{n}\cdot{\bm{q}}_{n})\right\|_{2}\left\|{\bm{V}}_{n}\right\|_{op}.$
(3)

This error bound is in line with the spectral errors studied in previous works
[27, 10].

We begin by observing that
$\mathrm{Attn}({\bm{q}}_{n},{\bm{K}}_{n},{\bm{V}}_{n})$ in ?? is the product
of the softmax vector $\mathtt{softmax}({\bm{K}}_{n}\cdot{\bm{q}}_{n})$ and
value matrix ${\bm{V}}_{n}$. This matrix-vector product can be approximated by
sub-sampling only $O(\varepsilon^{-2}d\log n)$ key-value pairs according to
the vector and matrix according to the squared norms of value tokens. This can
be implemented in a streaming setting using some variants of reservoir
sampling.

The other major computational challenge is computing the partition function in
the denominator of the softmax function, i.e.,
$\sum_{i\in[n]}\exp(\langle{\bm{k}}_{i},{\bm{q}}_{n}\rangle)$. To solve this,
we construct a data structure that can be stored in sublinear memory and
efficiently approximate
$\sum_{i\in[n]}\exp(\langle{\bm{k}}_{i},{\bm{q}}_{n}\rangle)$ up to
$1\pm\varepsilon$ factor for any query ${\bm{q}}_{n}$. Our method assumes that
the key tokens can be covered by a sublinear number of bounded diameter
clusters. This assumption is indeed weaker than the one made in [10], which in
the decoding setting translates to having key tokens belong to only one
cluster with a bounded diameter, while our approach allows for any sublinear
number of clusters. So, if the keys are composed of bounded diameter clusters
then we only need a small number of uniformly sampled keys from each cluster
to approximate the softmax normalizer efficiently and accurately. The central
task is to find these clusters in a streaming setting, and we achieve this
using an algorithm that is inspired by the streaming k-center algorithm of
[3].

In ?? and ?? we demonstrate that if the keys can be clustered into some
sublinear number $m=n^{1-\Omega(1)}$ of clusters with some bounded diameters,
then SubGen operates with sublinear
$O\left(\varepsilon^{-2}md\right)=O\left(\varepsilon^{-2}dn^{1-\Omega(1)}\right)$
memory and runtime and its output satisfies the approximation guarantee in ??.

In ??, we empirically compare SubGen to other KV cache compression methods
including the attention-score-based algorithm of [28] and the deterministic
eviction policy from [26]. Our results confirm that SubGen outperforms these
methods, particularly in question-answering tasks with various sequence
lengths.

##  2 Sublinear Time and Memory Algorithm

Our goal is to approximate the attention output in ?? with a space complexity
that is sublinear in context length $n$. To achieve this objective, we aim to
design the following data structure (DS) for efficiently approximating the
streaming attention mechanism:

###  2.1 Streaming Attention Data Structure

For every positive integer $n$ and every stream of token triplets
$({\bm{q}}_{1},{\bm{k}}_{1},{\bm{v}}_{1}),({\bm{q}}_{2},{\bm{k}}_{2},{\bm{v}}_{2}),\ldots({\bm{q}}_{n},{\bm{k}}_{n},{\bm{v}}_{n})$
where ${\bm{q}}_{i},{\bm{k}}_{i},{\bm{v}}_{i}\in\mathbb{R}^{d}$, we aim to
construct an efficient DS with the following properties:

  * •

The required memory space is sublinear in $n$, i.e., $o(n)$.

  * •

Upon the arrival of a new triplet
$({\bm{q}}_{n+1},{\bm{k}}_{n+1},{\bm{v}}_{n+1})$ in the stream, the time
complexity to update is sublinear in $n$, i.e., $o(n)$.

  * •

Given such data structure, there exists an algorithm that outputs an estimator
${\bm{z}}_{n}\in\mathbb{R}^{d}$ in sublinear time $o(n)$ such that:

$\displaystyle\left\|{\bm{z}}_{n}-\mathtt{softmax}({\bm{K}}_{n}\cdot{\bm{q}}_{n})^{\top}\cdot{\bm{V}}_{n}\right\|_{2}\leq\varepsilon\left\|\mathtt{softmax}({\bm{K}}_{n}\cdot{\bm{q}}_{n})\right\|_{2}\left\|{\bm{V}}_{n}\right\|_{op}.$
(4)

In the rest of this section, our focus is on developing an algorithm to
satisfy the above properties. Note that the attention output in ??, using the
definition of softmax, is equivalent to the following expression:

$\mathrm{Attn}({\bm{q}}_{n},{\bm{K}}_{n},{\bm{V}}_{n})=\frac{\exp({\bm{K}}_{n}\cdot{\bm{q}}_{n})^{\top}\cdot{\bm{V}}_{n}}{\sum_{i\in[n]}\exp(\langle{\bm{k}}_{i},{\bm{q}}_{n}\rangle)}.$

Thus, to compute the attention output we need to calculate:

  1. 1.

The matrix-vector product between ${\bm{V}}_{n}$ and
$\exp({\bm{K}}_{n}\cdot{\bm{q}}_{n})$.

  2. 2.

The partition function
$\sum_{i\in[n]}\exp(\langle{\bm{k}}_{i},{\bm{q}}_{n}\rangle)$.

Thus, our DS needs to efficiently approximate each of these two operations.
The matrix-vector product
$\exp({\bm{K}}_{n}\cdot{\bm{q}}_{n})^{\top}\cdot{\bm{V}}_{n}$ can be
approximated efficiently using standard sampling-based techniques.
Specifically, we make use of the row norm sampling approach [7, 6]. When
multiplying two matrices ${\bm{A}}\in\mathbb{R}^{m\times n}$ and
${\bm{B}}\in\mathbb{R}^{n\times p}$, we randomly sample an i.i.d. index
$i\in[n]$ with probability proportional to the $\ell_{2}$ norm of the $i$-th
row in ${\bm{B}}$. Then, we estimate ${\bm{A}}\cdot{\bm{B}}$ by the average of
the product between $i$-th column in ${\bm{A}}$ and $i$-th row in ${\bm{B}}$.
With this approximation, we need only $O(\varepsilon^{-2}d\log n)$ samples to
guarantee an $\varepsilon$ multiplicative error in spectral norm for
$\exp({\bm{K}}_{n}\cdot{\bm{q}}_{n})^{\top}\cdot{\bm{V}}_{n}$. Luckily, it can
be implemented in a streaming setting through a variant of reservoir sampling
[25].

The more challenging task is the sublinear-time approximation of the partition
function $\sum_{i\in[n]}\exp(\langle{\bm{k}}_{i},{\bm{q}}_{n}\rangle)$. We
construct a DS for computing this under the assumption that the keys in the
token stream are organized into a sublinear ($o(n)$) number of clusters. To be
more precise, we introduce the following notion of clusterability:

######  Definition 1 (Clusterability).

For a positive integer $m$ and a real-valued $\delta>0$, a dataset of points
${\bm{x}}_{1},{\bm{x}}_{2},\ldots{\bm{x}}_{n}\in\mathbb{R}^{d}$ is considered
$(m,\delta)$-clusterable if there exists a size-$m$ partition
$\mathcal{C}_{1},\mathcal{C}_{2},\ldots\mathcal{C}_{m}\subseteq\\{{\bm{x}}_{i}\\}_{i=1}^{n}$
of the dataset satisfying the following conditions:

  * •

$\mathcal{C}_{i}\cap\mathcal{C}_{j}=\emptyset$ for every $i\neq j$ and
$\bigcup_{j=1}^{m}\mathcal{C}_{j}=\\{{\bm{x}}_{i}\\}_{i=1}^{n}$.

  * •

for every $j\in[m]$ and every distinct pair
${\bm{y}},{\bm{z}}\in\mathcal{C}_{j}$,
$\left\|{\bm{y}}-{\bm{z}}\right\|_{2}\leq\delta$.

We demonstrate that under the assumption that the stream of keys
${\bm{k}}_{1},{\bm{k}}_{2},\ldots{\bm{k}}_{n}$ is $(m,\delta)$-clusterable as
defined in ??, with the number of clusters scaling sublinearly in stream
length ($m=o(n)$), it is possible to construct a DS with sublinear memory
space. The procedure for this DS is presented in ?? which we refer to as
SubGen.

To verify this in the practical settings, we plot key embeddings from open-
source LLMs in ?? and observe that they are indeed well clusterable on their
embedding space. This motivates us to utilize an efficient stream clustering
algorithm on key embeddings. In the remainder of this section, we provide a
detailed explanation for the execution of the algorithm while simultaneously
analyzing it through a series of lemmas.

1: inputs: stream of tokens $({\bm{q}}_{n},{\bm{k}}_{n},{\bm{v}}_{n})$ for
$n\in\mathbb{N}$, parameter $\delta>0$, positive integers $s,t$

2: initialize $\mu\leftarrow 0$, $\mathcal{D}\leftarrow\emptyset$,
$\mathcal{M}\leftarrow\begin{bmatrix}{\tt null},\stackrel{{\scriptstyle\times
s}}{{\cdots\cdots}}\end{bmatrix}$

3: repeat

4:
$\mathcal{D}\leftarrow\textsc{UpdateSoftmaxNormalizer}(\mathcal{D},\delta,t,{\bm{k}}_{n})$

5:
$\mathcal{M}\leftarrow\textsc{UpdateMatrixProduct}(\mathcal{M},s,\mu,{\bm{k}}_{n},{\bm{v}}_{n})$

6: $\mu\leftarrow\mu+\left\|{\bm{v}}_{n}\right\|_{2}^{2}$

7:
${\bm{z}}_{n}\leftarrow\textsc{QueryStreamAttn}(\mathcal{D},\mathcal{M},s,t,\mu,{\bm{q}}_{n})$

8: $n\leftarrow n+1$

9: output ${\bm{z}}_{n}$

10: until Token stream ends  Procedure UpdateSoftmaxNormalizer
($\mathcal{D},\delta,t,{\bm{k}}$)

11: Suppose input set
$\mathcal{D}=\\{({\bm{x}}_{i},\mathcal{S}_{i},n_{i}):i\in[m]\\}$

12: $i^{*}\leftarrow\mathop{\rm
arg~{}min}_{i\in[m]}\left\|{\bm{x}}_{i}-{\bm{k}}\right\|_{2}$

13: if $\left\|{\bm{k}}-{\bm{x}}_{i^{*}}\right\|_{2}\leq\delta$ then

14: $n_{i^{*}}\leftarrow n_{i^{*}}+1$

15: Suppose $\mathcal{S}_{i^{*}}$ is a list of $t$ vectors in $\mathbb{R}^{d}$

16: for $j\in[t]$ do

17: Flip a coin and with probability $p=\frac{1}{n_{i^{*}}}$, update the
$j^{th}$ entry of $\mathcal{S}_{i^{*}}$ as
$\mathcal{S}_{i^{*}}(j)\leftarrow{\bm{k}}$

18: end for

19: else

20:
$\mathcal{S}^{\prime}\leftarrow\begin{bmatrix}{\bm{k}},\stackrel{{\scriptstyle\times
t}}{{\cdots\cdots}}\end{bmatrix}$ (contains $t$ copies of ${\bm{k}}$)

21: $\mathcal{D}=\mathcal{D}\cup\\{({\bm{k}},\mathcal{S}^{\prime},1)\\}$

22: end if

123: return $\mathcal{D}$

2

Procedure UpdateMatrixProduct ($\mathcal{M},s,\mu,{\bm{k}},{\bm{v}}$)

24: Suppose $\mathcal{M}$ is a list of $s$ tuples of vectors in
$\mathbb{R}^{d}$

25: for $i\in[s]$ do

26: Flip a coin and with probability
$p=\frac{\left\|{\bm{v}}\right\|_{2}^{2}}{\mu+\left\|{\bm{v}}\right\|_{2}^{2}}$,
update the $i^{th}$ entry of $\mathcal{M}$ as
$\mathcal{M}(i)\leftarrow({\bm{k}},{\bm{v}})$

27: end for

328: return $\mathcal{M}$

4

Procedure QueryStreamAttn ($\mathcal{D},\mathcal{M},s,t,\mu,{\bm{q}}$)

29:
${\bm{z}}\leftarrow\sum_{({\bm{k}},{\bm{v}})\in\mathcal{M}}\frac{\mu}{s\cdot\left\|{\bm{v}}\right\|_{2}^{2}}\cdot\exp(\langle{\bm{q}},{\bm{k}}\rangle)\cdot{\bm{v}}$

30:
$\tau\leftarrow\sum_{({\bm{x}},\mathcal{S},n^{\prime})\in\mathcal{D}}\frac{n^{\prime}}{t}\cdot\sum_{{\bm{k}}\in\mathcal{S}}\exp(\langle{\bm{q}},{\bm{k}}\rangle)$

31: return ${\bm{z}}/\tau$

Algorithm 1 SubGen: Sublinear Streaming Attention

###  2.2 Matrix Product Data Structure

Here, we focus on the UpdateMatrixProduct primitive and establish its
correctness by introducing invariants that are maintained throughout the
stream processing. This primitive maintains and updates a list of $s$ elements
denoted by $\mathcal{M}$ in SubGen (??). Initially, this list is filled with
${\tt null}$ values. After processing the first token tuple
$({\bm{q}}_{1},{\bm{k}}_{1},{\bm{v}}_{1})$, this list is populated with $s$
copies of the first key and value $({\bm{k}}_{1},{\bm{v}}_{1})$. The procedure
UpdateMatrixProduct performs a variant of reservoir sampling upon observing
any new token in the stream. At any iteration $n$ of the stream, $\mathcal{M}$
is ensured to contain $s$ i.i.d. samples chosen at random from
$({\bm{k}}_{1},{\bm{v}}_{1}),({\bm{k}}_{2},{\bm{v}}_{2}),\ldots({\bm{k}}_{n},{\bm{v}}_{n})$
with probabilities proportional to $\left\|{\bm{k}}_{i}\right\|_{2}^{2}$. More
precisely, the following invariants hold:

######  Lemma 1 (Correctness of UpdateMatrixProduct).

For any positive integer $s$, at any iteration $n$ of the stream in 1 the
following properties are maintained:

  * •

$\mu=\sum_{i\in[n]}\left\|{\bm{v}}_{i}\right\|_{2}^{2}$.

  * •

$\mathcal{M}$ is a list of $s$ i.i.d. samples from
$\\{({\bm{k}}_{1},{\bm{v}}_{1}),({\bm{k}}_{2},{\bm{v}}_{2}),\ldots({\bm{k}}_{n},{\bm{v}}_{n})\\}$
where the probability distribution for each element $j\in[s]$ is
$\Pr[\mathcal{M}(j)=({\bm{k}}_{i},{\bm{v}}_{i})]=\frac{\left\|{\bm{v}}_{i}\right\|_{2}^{2}}{\sum_{l\in[n]}\left\|{\bm{v}}_{l}\right\|_{2}^{2}}$
for $i\in[n]$.

###### Proof.

The first property is trivial because $\mu$ is initialized at zero and is
updated in line 6 of the algorithm by adding the squared norms of
${\bm{v}}_{i}$’s. The proof of the second invariance is by induction. The base
of induction holds for $n=1$ because after processing the first token by
procedure UpdateMatrixProduct we have
$\Pr[\mathcal{M}(j)=({\bm{k}}_{1},{\bm{v}}_{1})]=\frac{\left\|{\bm{v}}_{1}\right\|_{2}^{2}}{\left\|{\bm{v}}_{1}\right\|_{2}^{2}}=1$
for $j\in[s]$.

Now suppose that the inductive hypothesis holds for $n$ and we prove it must
also hold for $n+1$. For any $j\in[s]$ in line 24 of ?? with probability
$p=\frac{\left\|{\bm{v}}_{n+1}\right\|_{2}^{2}}{\mu+\left\|{\bm{v}}_{n+1}\right\|_{2}^{2}}$,
$\mathcal{M}(j)$ gets updated to $({\bm{k}}_{n+1},{\bm{v}}_{n+1})$. Since we
showed that $\mu=\sum_{i\in[n]}\left\|{\bm{v}}_{i}\right\|_{2}^{2}$ we have:

$\Pr[\mathcal{M}(j)=({\bm{k}}_{n+1},{\bm{v}}_{n+1})]=\frac{\left\|{\bm{v}}_{n+1}\right\|_{2}^{2}}{\sum_{l\in[n+1]}\left\|{\bm{v}}_{l}\right\|_{2}^{2}}.$

Moreover with probability
$1-p=\frac{\mu}{\mu+\left\|{\bm{v}}_{n+1}\right\|_{2}^{2}}$, $\mathcal{M}(j)$
keeps its previous value. Using the inductive hypothesis we have that for
every $i\in[n]$:

$\displaystyle\Pr[\mathcal{M}(j)=({\bm{k}}_{i},{\bm{v}}_{i})]$
$\displaystyle=\frac{\left\|{\bm{v}}_{i}\right\|_{2}^{2}}{\sum_{l\in[n]}\left\|{\bm{v}}_{l}\right\|_{2}^{2}}\cdot\frac{\sum_{l\in[n]}\left\|{\bm{v}}_{l}\right\|_{2}^{2}}{\sum_{l\in[n+1]}\left\|{\bm{v}}_{l}\right\|_{2}^{2}}=\frac{\left\|{\bm{v}}_{i}\right\|_{2}^{2}}{\sum_{l\in[n+1]}\left\|{\bm{v}}_{l}\right\|_{2}^{2}}.$

This completes the proof. ∎

###  2.3 Softmax Normalizer (Partition Function) DS

Here we delve into a detailed discussion of the UpdateSoftmaxNormalizer
primitive. This primitive constructs and maintains a DS denoted by
$\mathcal{D}$, enabling accurate approximation of the partition function in
the softmax denominator for any query. A crucial requirement for the
efficiency of this primitive is that the key tokens must be
$(m,\delta)$-clusterable, as per ??. Our algorithm locates and stores a
subsampled representation of each cluster in $\mathcal{D}$ in a small memory.
Particularly, to achieve sublinear memory complexity, instead of keeping all
keys in each cluster which would require $O(n)$ memory space, we maintain only
a random subset of $t$ samples from each cluster.

Initially, $\mathcal{D}$ is an empty set. As new tokens in the stream are
processed, new clusters get added to this set. Each cluster is characterized
by a representative point, which is the first key assigned to that cluster by
our algorithm. Throughout stream processing, we compute the distance between
the new key token and each existing cluster. Here the distance to an existing
cluster is defined as the distance to the aforementioned representative of the
cluster. If there is a cluster whose distance is less than $\delta$, then the
token is assigned to the nearest cluster, and we update our random samples of
keys from this cluster using reservoir sampling. If the distance from all
existing clusters is more than $\delta$, we introduce a new cluster in
$\mathcal{D}$, and the new key becomes the representative of this new cluster.
At any point in the stream, this algorithm identifies at most $m$ clusters if
the keys so far are $(m,\delta)$-clusterable. If $m$ grows sublinearly in the
stream length $n$, the memory and update time of our algorithm will be
sublinear as well. Formally, we prove that the following invariant holds:

######  Lemma 2 (Correctness of UpdateSoftmaxNormalizer).

For any $\delta>0$, any positive integer $t$, at any iteration $n$ of the
stream in ?? the following properties are maintained. $\mathcal{D}$ is a set
of $m$ items of the form
$\mathcal{D}=\left\\{({\bm{x}}_{i},\mathcal{S}_{i},n_{i}):i\in[m]\right\\}$,
where there exists a partition of keys into $m$ disjoint subsets
$\mathcal{C}_{1},\mathcal{C}_{2},\ldots\mathcal{C}_{m}\subseteq\\{{\bm{k}}_{i}\\}_{i=1}^{n}$
satisfying $\bigcup_{j=1}^{m}\mathcal{C}_{j}=\\{{\bm{k}}_{i}\\}_{i=1}^{n}$ and
$\mathcal{C}_{i}\cap\mathcal{C}_{j}=\emptyset$ for every $i\neq j$, such that
for every $i\in[m]$:

  1. 1.

${\bm{x}}_{i}\in\mathcal{C}_{i}$,

  2. 2.

$n_{i}=|\mathcal{C}_{i}|$,

  3. 3.

$\left\|{\bm{x}}_{i}-{\bm{k}}^{\prime}\right\|_{2}\leq\delta$ for every
${\bm{k}}^{\prime}\in\mathcal{C}_{i}$,

  4. 4.

$\left\|{\bm{x}}_{i}-{\bm{x}}_{j}\right\|_{2}>\delta$ for every $i\neq j$,

  5. 5.

$\mathcal{S}_{i}$ is a set of $t$ i.i.d. uniform samples from the set
$\mathcal{C}_{i}$.

###### Proof.

The proof is by induction on the stream length $n$. The base of induction
trivially holds for $n=0$, where $\mathcal{D}$ is an empty set. To prove the
inductive step suppose that the inductive hypothesis holds for some $n$.
Specifically, suppose that $\mathcal{D}$ is a set of $m$ items of the form
$\mathcal{D}=\left\\{({\bm{x}}_{i},\mathcal{S}_{i},n_{i}):i\in[m]\right\\}$
and there exists a partition of keys into $m$ disjoint subsets
$\mathcal{C}_{1},\mathcal{C}_{2},\ldots\mathcal{C}_{m}\subseteq\\{{\bm{k}}_{i}\\}_{i=1}^{n}$
as per in the lemma statement, such that for every $i\in[m]$: (1)
${\bm{x}}_{i}\in\mathcal{C}_{i}$, (2) $n_{i}=|\mathcal{C}_{i}|$, (3)
$\left\|{\bm{x}}_{i}-{\bm{k}}^{\prime}\right\|_{2}\leq\delta$ for every
${\bm{k}}^{\prime}\in\mathcal{C}_{i}$, (4)
$\left\|{\bm{x}}_{i}-{\bm{x}}_{j}\right\|_{2}>\delta$ for every $i\neq j$, and
(5) $\mathcal{S}_{i}$ is a set of $t$ i.i.d. uniform samples from the set
$\mathcal{C}_{i}$. Given this assumption, we prove that the inductive step
also holds for after processing the $(n+1)$-th key in the stream
${\bm{k}}_{n+1}$.

In the next iteration, specifically in line 12 of UpdateSoftmaxNormalizer, the
algorithm finds the index $i^{*}\in[m]$ such that
$\left\|{\bm{x}}_{i^{*}}-{\bm{k}}_{n+1}\right\|_{2}$ is minimized. Two cases
arise:

#### Case 1: $\left\|{\bm{x}}_{i^{*}}-{\bm{k}}_{n+1}\right\|_{2}\leq\delta$.

In this case, the algorithm increments $n_{i^{*}}\leftarrow n_{i^{*}}+1$ in
line 14. Consider the new partitioning of the keys defined as
$\mathcal{C}_{i}^{\prime}=\mathcal{C}_{i}$ for $i\neq i^{*}$ and
$\mathcal{C}_{i^{*}}^{\prime}=\mathcal{C}_{i^{*}}\cup\\{{\bm{k}}_{n+1}\\}$. It
follows from the inductive hypothesis that for every $i\in[m]$: (1)
${\bm{x}}_{i}\in\mathcal{C}_{i}^{\prime}$, (2)
$n_{i}=|\mathcal{C}_{i}^{\prime}|$, (3)
$\left\|{\bm{x}}_{i}-{\bm{k}}^{\prime}\right\|_{2}\leq\delta$ for every
${\bm{k}}^{\prime}\in\mathcal{C}_{i}^{\prime}$, and (4)
$\left\|{\bm{x}}_{i}-{\bm{x}}_{j}\right\|_{2}>\delta$ for every $i\neq j$ hold
after the $n+1$-th iteration. Furthermore, since the algorithm does not alter
the lists $\mathcal{S}_{i}$ for $i\neq i^{*}$, we have that (5)
$\mathcal{S}_{i}$ is a set of $t$ i.i.d. uniform samples from the set
$\mathcal{C}_{i}^{\prime}$ for any $i\neq i^{*}$. On the other hand, the
algorithm in line 17 performs reservoir sampling on the set
$\mathcal{S}_{i^{*}}$ with new element ${\bm{k}}_{n+1}$ which implies that
$\mathcal{S}_{i^{*}}$ is a set of $t$ i.i.d. uniform samples from the set
$\mathcal{C}_{i^{*}}^{\prime}$. This completes the inductive step in the first
case.

#### Case 2: $\left\|{\bm{x}}_{i^{*}}-{\bm{k}}_{n+1}\right\|_{2}>\delta$.

In this case, the algorithm adds a new element to $\mathcal{D}$, thus, the
updated set is
$\mathcal{D}^{\prime}=\\{({\bm{x}}_{i},\mathcal{S}_{i},n_{i}):i\in[m+1]\\}$
with ${\bm{x}}_{m+1}={\bm{k}}_{n+1}$ and $n_{m+1}=1$. If we consider the new
partitioning of keys to be
$\mathcal{C}_{1},\mathcal{C}_{2},\ldots\mathcal{C}_{m},\mathcal{C}_{m+1}$,
where $\mathcal{C}_{m+1}=\\{{\bm{k}}_{n+1}\\}$, we can use the inductive
hypothesis to deduce that for any $i\in[m+1]$: (1)
${\bm{x}}_{i}\in\mathcal{C}_{i}$, (2) $n_{i}=|\mathcal{C}_{i}|$, (3)
$\left\|{\bm{x}}_{i}-{\bm{k}}^{\prime}\right\|_{2}\leq\delta$ for every
${\bm{k}}^{\prime}\in\mathcal{C}_{i}$, and (4)
$\left\|{\bm{x}}_{i}-{\bm{x}}_{j}\right\|_{2}>\delta$ for every $i\neq j$ hold
after the $n+1$-th iteration of the stream. Furthermore, $\mathcal{S}_{m+1}$
is defined to be a list of $t$ copies of ${\bm{k}}_{n+1}$, thus, (5)
$\mathcal{S}_{i}$ is a set of $t$ i.i.d. uniform samples from the set
$\mathcal{C}_{i}$ for any $i\in[m+1]$. This completes the inductive step in
this case and also concludes the proof of the lemma. ∎

###  2.4 Streaming Attention: Main Theorem

Now we are ready to analyze the end-to-end performance of SubGen and prove the
main theorem. We show that, given the data structures created throughout the
stream and analyzed in ?? and ??, the primitive QueryStreamAttn can
efficiently output an accurate approximation to the streaming attention,
satisfying ??.

Our analysis unfolds in two steps. First, we establish that the data
structures created by UpdateSoftmaxNormalizer and UpdateMatrixProduct can be
stored in small memory and updated very quickly if the sequence of keys is
clusterable into a sublinear number of clusters. Then we show that the
QueryStreamAttn can use these data structures to produce an accurate attention
output for any given query. Our main result is as follows:

######  Theorem 1 (Efficiency and Correctness of ??).

For any $\delta,r,\varepsilon>0$, any positive integers $n,d$, and any
sequence of tokens
$({\bm{q}}_{1},{\bm{k}}_{1},{\bm{v}}_{1}),({\bm{q}}_{2},{\bm{k}}_{2},{\bm{v}}_{2}),\ldots({\bm{q}}_{n},{\bm{k}}_{n},{\bm{v}}_{n})$
where ${\bm{q}}_{i},{\bm{k}}_{i},{\bm{v}}_{i}\in\mathbb{R}^{d}$, suppose that
the followings hold

  * •

$t=\Omega\left(\varepsilon^{-2}\cdot e^{2\delta\cdot r}\log n\right)$,

  * •

$s=\Omega(\varepsilon^{-2}\cdot d)$,

  * •

$\left\|{\bm{q}}_{n}\right\|_{2}\leq r$.

Then, SubGen (??) at $n$-th step of the stream processing outputs a vector
${\bm{z}}_{n}\in\mathbb{R}^{d}$ that satisfies ?? with probability at least
$0.99$. Furthermore, if the keys
${\bm{k}}_{1},{\bm{k}}_{2},\ldots{\bm{k}}_{n}$ are $(m,\delta)$-clusterable as
per ??, then both the total memory of the algorithm and its runtime during the
$n$-th iteration is bounded by $O(d\cdot(mt+s))$.

###### Proof.

We start the correctness proof by observing that all preconditions of ?? are
satisfied, allowing us to invoke this lemma. Let the partition of keys into
disjoint subsets be denoted by
$\mathcal{C}_{1},\mathcal{C}_{2},\ldots\mathcal{C}_{m^{\prime}}\subseteq\\{{\bm{k}}_{i}\\}_{i=1}^{n}$
satisfying
$\bigcup_{j=1}^{m^{\prime}}\mathcal{C}_{j}=\\{{\bm{k}}_{i}\\}_{i=1}^{n}$ and
$\mathcal{C}_{i}\cap\mathcal{C}_{j}=\emptyset$ for every $i\neq j$ as per ??
for some positive integer $m^{\prime}$. Rewriting the partition function in
the attention denominator gives:

$\sum_{j\in[n]}\exp(\langle{\bm{k}}_{j},{\bm{q}}_{n}\rangle)=\sum_{i\in[m^{\prime}]}\sum_{{\bm{k}}^{\prime}\in\mathcal{C}_{i}}\exp(\langle{\bm{k}}^{\prime},{\bm{q}}_{n}\rangle).$

Now by property (3) in ?? and triangle inequality, for every
$i\in[m^{\prime}]$ and every
${\bm{k}}^{\prime},{\bm{k}}^{\prime\prime}\in\mathcal{C}_{i}$ we have:

$\left\|{\bm{k}}^{\prime}-{\bm{k}}^{\prime\prime}\right\|_{2}\leq\left\|{\bm{k}}^{\prime}-{\bm{x}}_{i}\right\|_{2}+\left\|{\bm{k}}^{\prime\prime}-{\bm{x}}_{i}\right\|_{2}\leq
2\delta.$

Therefore, using the precondition of the theorem on
$\left\|{\bm{q}}_{n}\right\|_{2}\leq r$ we have

$\exp(\langle{\bm{k}}^{\prime},{\bm{q}}_{n}\rangle)/\exp(\langle{\bm{k}}^{\prime\prime},{\bm{q}}_{n}\rangle)\leq
e^{2\delta\cdot r}.$

Using the above inequality and the assumption in the theorem statement
regarding $t=\Omega\left(\varepsilon^{-2}\cdot e^{2\delta\cdot r}\log
n\right)$ combined with the properties (2) and (5) proved in ??, we can invoke
Chernoff-Hoeffding inequality (see e.g., [16]) along with union bound to
conclude that the following holds simultaneously for all $i\in[m^{\prime}]$
with probability at least $1-\frac{1}{\text{poly}(n)}$:

$\frac{n_{i}}{t}\cdot\sum_{{\bm{k}}^{\prime}\in\mathcal{S}_{i}}\exp(\langle{\bm{q}}_{n},{\bm{k}}^{\prime}\rangle)\in(1\pm\varepsilon/3)\cdot\sum_{{\bm{k}}^{\prime}\in\mathcal{C}_{i}}\exp(\langle{\bm{k}}^{\prime},{\bm{q}}_{n}\rangle)$

Since the terms above are positive, by summing up the given inequality for all
$i\in[m^{\prime}]$, we find that the quantity $\tau$ computed in line 27 of ??
satisfies the following:

$\Pr\left[\tau\in(1\pm\varepsilon/3)\sum_{j\in[n]}\exp(\langle{\bm{k}}_{j},{\bm{q}}_{n}\rangle)\right]\geq
0.995$ (5)

Next, we invoke ?? to derive an error bound on the approximate matrix-vector
product between the softmax vector and the matrix of values ${\bm{V}}_{n}$. By
leveraging well-established techniques in approximate matrix products, such as
the standard result from [7], and using the conclusion of ?? regarding
$\mathcal{M}$ as a list of $s=\Omega(\varepsilon^{-2}\cdot d)$ i.i.d. sample
from the probability distribution
$\Pr[\mathcal{M}(j)=({\bm{k}}_{i},{\bm{v}}_{i})]=\frac{\left\|{\bm{v}}_{i}\right\|_{2}^{2}}{\sum_{l\in[n]}\left\|{\bm{v}}_{l}\right\|_{2}^{2}}$
for $i\in[n]$ for $i\in[n]$ and $j\in[s]$, we have that vector ${\bm{z}}$
computed in line 26 of ?? satisfies the following inequality with a
probability of at least $0.995$:

$\displaystyle\left\|{\bm{z}}-\exp({\bm{K}}_{n}\cdot{\bm{q}}_{n})^{\top}\cdot{\bm{V}}_{n}\right\|_{2}\leq\frac{\varepsilon}{3}\left\|\exp({\bm{K}}_{n}\cdot{\bm{q}}_{n})\right\|_{2}\|{\bm{V}}_{n}\|_{op}$
(6)

Now by combining inequalities in ?? and ?? using union bound and triangle
inequality we find that the output of ?? computed in line 28 as
${\bm{z}}/\tau$ satisfies the following with probability at least $0.99$

$\displaystyle\left\|{\bm{z}}/\tau-\mathtt{softmax}({\bm{K}}_{n}\cdot{\bm{q}}_{n})^{\top}\cdot{\bm{V}}_{n}\right\|_{2}\leq\varepsilon\left\|\mathtt{softmax}({\bm{K}}_{n}\cdot{\bm{q}}_{n})\right\|_{2}\|{\bm{V}}_{n}\|_{op}.$

This completes the correctness proof.

Figure 1: A t-SNE plot of cached keys (first row) and values (second row)
embeddings over $1024$ timesteps from Llama2-7B using MT Bench dataset. We
pick $\ell$-layer where $\ell\in\\{0,7,15,23,31\\}$ and head IDs are chosen
uniformly at random. Key embeddings are more clusterable than value ones. The
green dots represent the centers from the greedy k-center algorithm [8] where
k=$16$.

#### Memory and Runtime.

First, note that the memory requirement for storing the list $\mathcal{M}$ in
?? is $O(sd)$ because it contains $s$ pairs of $d$-dimensional vectors. Next,
to bound the memory requirement for storing $\mathcal{D}$ we need to bound the
size of this set which we denoted by $m^{\prime}$. According to properties (1)
and (4) in ??, for every $i\in[m^{\prime}]$ there exist
${\bm{x}}_{i}\in\\{{\bm{k}}_{1},{\bm{k}}_{2},\ldots{\bm{k}}_{n}\\}$ such that
$\left\|{\bm{x}}_{i}-{\bm{x}}_{j}\right\|_{2}>\delta$ for $i\neq j$. Given the
assumption in the theorem statement that keys are $(m,\delta)$-clusterable, by
the definition of clusterability in ?? along with the pigeonhole principle, we
must have $m^{\prime}\leq m$. Therefore storing $\mathcal{D}$ will require
$O(m^{\prime}td)=O(mtd)$ because it is a set of $m^{\prime}$ elements, and
each element of this set is a list of $t$ vectors in dimension $d$.

Three major operations dominate the runtime of the $n$-th iteration. Firstly,
executing UpdateSoftmaxNormalizer requires computing $m^{\prime}$ distances in
line 12 that takes $O(md)$ time. Additionally, the for loop in line 16 takes
$O(td)$ time. Secondly, UpdateMatrixProduct has a runtime bounded by $O(sd)$.
Thirdly, running QueryStreamAttn involves $O(sd)$ operations in line 26 and
$O(m^{\prime}td)=O(mtd)$ operations in line 27. As a result, the total runtime
of ?? in $n$-th iteration is $O(mtd+sd)$. ∎

?? demonstrates that if the keys can be clustered into some sublinear number
$m=n^{1-\Omega(1)}$ of clusters with diameters at most $\delta$, and the
queries have bounded $\ell_{2}$ norms of at most $r$ such that the product of
the cluster diameter and maximum $\ell_{2}$ norm of queries is bounded by
$\delta r=o(\log n)$, then ?? operates with sublinear
$O\left(\varepsilon^{-2}\cdot mdn^{o(1)}\right)=O\left(\varepsilon^{-2}\cdot
dn^{1-\Omega(1)}\right)$ memory and runtime. We summarize this in the
following corollary:

######  Corollary 1.

Suppose the preconditions of ?? hold. If the diameter of key token clusters
$\delta$ and the maximum $\ell_{2}$ norm of queries $r$ satisfy $\delta
r=o(\log n)$, then the total memory and runtime of ?? are bounded by
$O\left(\varepsilon^{-2}\cdot dmn^{o(1)}\right)$. Moreover, if the number of
key token clusters $m$ grows as a sublinear function of $n$, i.e., as
$m=n^{1-\Omega(1)}$, then the memory and runtime are bounded by
$O\left(\varepsilon^{-2}\cdot dn^{1-\Omega(1)}\right)$.

|  $n=$ 5k |  $n=$ 7k |  $n=$ 9k
---|---|---|---
Algorithm | Cache Size (GB) | Accuracy | Cache Size (GB) | Accuracy | Cache Size (GB) | Accuracy
Exact | 2.351 | 0.98 | 3.488 | 1.0 | 4.613 | 0.68
Sink [26] | 1.511 (35% $\downarrow$) | 0.56 | 2.012 (42% $\downarrow$) | 0.56 | 2.262 (50% $\downarrow$) | 0.38
H2O [28] | 1.511 (35% $\downarrow$) | 0.66 | 2.012 (42% $\downarrow$) | 0.58 | 2.262 (50% $\downarrow$) | 0.38
SubGen (this work) | 1.512 (35% $\downarrow$) | 0.86 | 2.012 (42% $\downarrow$) | 0.66 | 2.262 (50% $\downarrow$) | 0.44

Table 1: Results on accuracy of line retrieval from LongEval [13] dataset with
context length 5k-9k. Under the sublinear budgets on cache size, the proposed
algorithm based on k-center algorithm outperforms other methods over all
sequence lengths.

##  3 Experiments

In this section, we report the empirical results of the proposed algorithm
with memory footprint reduction and performance on question-answering
benchmark datasets. For all experiments, we use a single NVIDIA RTX6000 GPU
with 48 GB VRAM.

###  3.1 Ablation Study on Clusterability

We first demonstrate that cached embeddings over long token generations are
indeed well clusterable. To this end, we collect key and value embeddings from
Llama-2-7B [23] with MT Bench dataset [29] while the model generates a
sequence of $1024$ tokens. We then visualize the embeddings using t-SNE [24]
across various layers and heads, identifying cluster center points through the
greedy k-center algorithm [8].

As illustrated in ??, our observations reveal that key embeddings (first row)
exhibit a higher degree of clusterability compared to value embeddings across
all randomly selected layers and heads. Furthermore, we note that the cluster
centers (indicated by green dots) corresponding to the key embeddings are
evenly distributed across the entire embedding space. In particular, the key
embeddings demonstrate significant dispersion across different time steps, and
their cluster centers are distributed over the entire embedding space.

This behavior comes from the use of Rotary Position Embedding (RoPE) [21] in
Llama-2-type models which introduces rotational transformations to both query
and key embeddings based on their relative positions. Hence, the key
embeddings appear to be well-separated in their projected space, while the
values show an unstructured and random distribution within their space. These
findings serve as a motivation for the development of an efficient key-value
(KV) compression algorithm that leverages the clustering properties of key
embeddings.

###  3.2 End-to-end Text Generation

We next evaluate our proposed algorithm on long-context line retrieval task in
LongEval [13]333https://github.com/DachengLi1/LongChat/blob/longeval
benchmark. The task involves long-context line retrieval from extensive
documents, each comprising multiple lines, complete with line numbers and
topics. The objective is to precisely retrieve a specified number of lines
corresponding to a target topic. We vary the number of lines, representing the
number of targets, to 200, 300, and 400 and they correspond to sequence
lengths of $n=$5,000, 7,000, and 9,000, respectively. Each dataset contains 50
distinct questions, and we systematically extract the number from the
generated answers and compute accuracies. The answers are generated employing
the longchat-7B model444https://huggingface.co/lmsys/longchat-7b-v1.5-32k,
which is a fine-tuned version of the Llama-2-7B model with long-range context
length.

We compare our method to two KV cache compression algorithms; H2O [28], which
retains cached tokens with high cumulative attention scores, and Attention
Sink [26], a method that deterministically selects some initial and recent
tokens. Specifically, both of these prior works have highlighted the
significance of recent token embeddings in generating meaningful responses. To
leverage this insight, we integrate it with our clustering approach. More
precisely, our strategy consistently retains the most recent $r$ token
embeddings, in addition to $k$ centers selected from the remaining tokens. In
a streaming context, this strategy is often referred to as a sliding window.
We apply the greedy k-center clustering algorithm once to compress the entire
KV caches. To make comparisons fair, we set cache memory budgets of all
algorithms identical (i.e., $r+k$), which scales sublinearly with the context
length denoted as $n$.

The results are reported in ??. We observe that our clustering-based method
consistently outperforms other algorithms across all sequence lengths. For
instance, we achieve an accuracy of 44% while utilizing only half of the
cached KV embeddgins with a length of 9k tokens, whereas both H2O and Sink can
achieve accuracies 10% lower. This finding suggests that maintaining the
embedding information holds greater significance in sustaining the performance
of LLMs compared to attention scores and positional information.

##  4 Conclusion

In this work, we develop SubGen, an efficient KV cache compression algorithm
via stream clustering. Our motivation is that cached keys are well clusterable
in their embedding space and we apply a greedy-type clustering algorithm to
find the most representative embeddings. Under assumptions on bounded query
norm and clusterability, we analyze that our algorithm can guarantee a
spectral error bound with sublinear time and memory. We further integrate
keeping recent tokens to the proposed clustering approach. For zero-shot line
retrieval tasks, our algorithm outperforms other KV cache compression
algorithms with the same memory budget.

## References

  * Achiam et al. [2023] Achiam, J., Adler, S., Agarwal, S., Ahmad, L., Akkaya, I., Aleman, F. L., Almeida, D., Altenschmidt, J., Altman, S., Anadkat, S., et al.  Gpt-4 technical report.  _arXiv preprint arXiv:2303.08774_ , 2023.
  * Alman & Song [2023] Alman, J. and Song, Z.  Fast attention requires bounded entries.  _arXiv preprint arXiv:2302.13214_ , 2023.
  * Charikar et al. [1997] Charikar, M., Chekuri, C., Feder, T., and Motwani, R.  Incremental clustering and dynamic information retrieval.  In _Symposium on the Theory of Computing (STOC)_ , 1997.
  * Chen et al. [2021] Chen, M., Tworek, J., Jun, H., Yuan, Q., Pinto, H. P. d. O., Kaplan, J., Edwards, H., Burda, Y., Joseph, N., Brockman, G., et al.  Evaluating large language models trained on code.  _arXiv preprint arXiv:2107.03374_ , 2021.
  * Chiang et al. [2023] Chiang, W.-L., Li, Z., Lin, Z., Sheng, Y., Wu, Z., Zhang, H., Zheng, L., Zhuang, S., Zhuang, Y., Gonzalez, J. E., et al.  Vicuna: An open-source chatbot impressing gpt-4 with 90%* chatgpt quality.  2023\.  URL https://vicuna.lmsys.org.
  * Cohen et al. [2016] Cohen, M. B., Nelson, J., and Woodruff, D. P.  Optimal Approximate Matrix Product in Terms of Stable Rank.  In _International Colloquium on Automata, Languages, and Programming (ICALP)_ , 2016.
  * Drineas & Kannan [2001] Drineas, P. and Kannan, R.  Fast Monte-Carlo algorithms for approximate matrix multiplication.  In _Foundations of Computer Science (FOCS)_ , 2001.
  * Dyer & Frieze [1985] Dyer, M. E. and Frieze, A. M.  A simple heuristic for the p-centre problem.  _Operations Research Letters_ , 1985.
  * Ge et al. [2023] Ge, S., Zhang, Y., Liu, L., Zhang, M., Han, J., and Gao, J.  Model Tells You What to Discard: Adaptive KV Cache Compression for LLMs.  _arXiv preprint arXiv:2310.01801_ , 2023.
  * Han et al. [2024] Han, I., Jayaram, R., Karbasi, A., Mirrokni, V., Woodruff, D. P., and Zandieh, A.  Hyperattention: Long-context attention in near-linear time.  _International Conference on Learning Representations (ICLR)_ , 2024\.
  * Ho et al. [2022] Ho, J., Chan, W., Saharia, C., Whang, J., Gao, R., Gritsenko, A., Kingma, D. P., Poole, B., Norouzi, M., Fleet, D. J., et al.  Imagen video: High definition video generation with diffusion models.  _arXiv preprint arXiv:2210.02303_ , 2022.
  * Hong et al. [2023] Hong, K., Dai, G., Xu, J., Mao, Q., Li, X., Liu, J., Chen, K., Dong, H., and Wang, Y.  Flashdecoding++: Faster large language model inference on gpus.  _arXiv preprint arXiv:2311.01282_ , 2023.
  * Li et al. [2023] Li, D., Shao, R., Xie, A., Sheng, Y., Zheng, L., Gonzalez, J., Stoica, I., Ma, X., and Zhang, H.  How Long Can Context Length of Open-Source LLMs truly Promise? In _NeurIPS 2023 Workshop on Instruction Tuning and Instruction Following_ , 2023.
  * Liu et al. [2023a] Liu, Z., Desai, A., Liao, F., Wang, W., Xie, V., Xu, Z., Kyrillidis, A., and Shrivastava, A.  Scissorhands: Exploiting the Persistence of Importance Hypothesis for LLM KV Cache Compression at Test Time.  _Neural Information Processing Systems (NeurIPS)_ , 2023a.
  * Liu et al. [2023b] Liu, Z., Wang, J., Dao, T., Zhou, T., Yuan, B., Song, Z., Shrivastava, A., Zhang, C., Tian, Y., Re, C., et al.  Deja vu: Contextual sparsity for efficient llms at inference time.  In _International Conference on Machine Learning (ICML)_ , 2023b.
  * McDiarmid [1998] McDiarmid, C.  Concentration.  In _Probabilistic methods for algorithmic discrete mathematics_. Springer, 1998.
  * Radford et al. [2021] Radford, A., Kim, J. W., Hallacy, C., Ramesh, A., Goh, G., Agarwal, S., Sastry, G., Askell, A., Mishkin, P., Clark, J., et al.  Learning transferable visual models from natural language supervision.  In _International Conference on Machine Learning (ICML)_ , 2021.
  * Roziere et al. [2023] Roziere, B., Gehring, J., Gloeckle, F., Sootla, S., Gat, I., Tan, X. E., Adi, Y., Liu, J., Remez, T., Rapin, J., et al.  Code llama: Open foundation models for code.  _arXiv preprint arXiv:2308.12950_ , 2023.
  * Sarlos et al. [2023] Sarlos, T., Song, X., Woodruff, D., and Zhang, Q.  Hardness of low rank approximation of entrywise transformed matrix products.  In _Thirty-seventh Conference on Neural Information Processing Systems_ , 2023.
  * Sheng et al. [2023] Sheng, Y., Zheng, L., Yuan, B., Li, Z., Ryabinin, M., Chen, B., Liang, P., Ré, C., Stoica, I., and Zhang, C.  Flexgen: High-throughput generative inference of large language models with a single gpu.  In _International Conference on Machine Learning (ICML)_ , 2023.
  * Su et al. [2024] Su, J., Ahmed, M., Lu, Y., Pan, S., Bo, W., and Liu, Y.  Roformer: Enhanced transformer with rotary position embedding.  _Neurocomputing_ , 2024.
  * Taori et al. [2023] Taori, R., Gulrajani, I., Zhang, T., Dubois, Y., Li, X., Guestrin, C., Liang, P., and Hashimoto, T. B.  Stanford alpaca: An instruction-following llama model.  https://github.com/tatsu-lab/stanford_alpaca, 2023.
  * Touvron et al. [2023] Touvron, H., Martin, L., Stone, K., Albert, P., Almahairi, A., Babaei, Y., Bashlykov, N., Batra, S., Bhargava, P., Bhosale, S., et al.  Llama 2: Open foundation and fine-tuned chat models.  _arXiv preprint arXiv:2307.09288_ , 2023.
  * Van der Maaten & Hinton [2008] Van der Maaten, L. and Hinton, G.  Visualizing data using t-SNE. _Journal of Machine Learning Research (JMLR)_ , 2008.
  * Vitter [1985] Vitter, J. S.  Random sampling with a reservoir.  _ACM Transactions on Mathematical Software (TOMS)_ , 1985.
  * Xiao et al. [2023] Xiao, G., Tian, Y., Chen, B., Han, S., and Lewis, M.  Efficient streaming language models with attention sinks.  _arXiv preprint arXiv:2309.17453_ , 2023.
  * Zandieh et al. [2023] Zandieh, A., Han, I., Daliri, M., and Karbasi, A.  Kdeformer: Accelerating transformers via kernel density estimation.  In _International Conference on Machine Learning (ICML)_ , 2023.
  * Zhang et al. [2023] Zhang, Z., Sheng, Y., Zhou, T., Chen, T., Zheng, L., Cai, R., Song, Z., Tian, Y., Ré, C., Barrett, C., et al.  H2O: Heavy-Hitter Oracle for Efficient Generative Inference of Large Language Models.  _Neural Information Processing Systems (NeurIPS)_ , 2023.
  * Zheng et al. [2023] Zheng, L., Chiang, W.-L., Sheng, Y., Zhuang, S., Wu, Z., Zhuang, Y., Lin, Z., Li, Z., Li, D., Xing, E., et al.  Judging LLM-as-a-judge with MT-Bench and Chatbot Arena.  _arXiv preprint arXiv:2306.05685_ , 2023.  URL https://github.com/lm-sys/FastChat/tree/main/fastchat/llm_judge.


# Virtual Multi-Topology Routing for QoS Constraints

Nicolas Huin1, Sébastien Martin2, Jérémie Leguay2
1IMT Atlantique, IRISA UMR CNRS 6074, F-35700 Rennes, France 2Huawei
Technologies Ltd., Paris Research Center, France.

###### Abstract

Multi-topology routing (MTR) provides an attractive alternative to segment
routing for traffic engineering when network devices cannot be upgraded.
However, due to a high overhead in terms of link state messages exchanged by
topologies and the need to frequently update link weights to follow evolving
network conditions, MTR is often limited to a small number of topologies and
the satisfaction of loose QoS constraints. To overcome these limitations we
propose vMTR, an MTR extension where demands are routed over virtual
topologies that are silent, i.e., they do not exchange LSA messages, and that
are continuously derived from a very limited set of real topologies,
optimizing each a QoS parameter. In this context, we present a polynomial and
exact algorithm for vMTR and, as a benchmark, a local search algorithm for
MTR. We show that vMTR helps reducing drastically the number of real
topologies and that it is more robust to QoS changes.

###### Index Terms:

Multi-Topology Routing, Traffic Engineering, Quality of Service, Segment
Routing.

##  I Introduction

Traffic Engineering (TE) policies [1] are key to controlling how bandwidth is
shared in IP (Internet Protocol) networks. These policies can be tuned to
satisfy application requirements in terms of end-to-end delay, packet loss,
and jitter, for instance. They can be optimized by a centralized controller
and applied by routers to steer traffic with data plane mechanisms such as
Segment Routing (SR) [2, 3] or MPLS (Multi-Protocol Label Switching). However,
most IP networks today still rely on shortest path and hop-by-hop routing from
Interior Gateway Protocols (IGP).

Extensions of IGP protocols such as Multi-Topology Routing (MTR) [4, 5] or
FlexAlgo [6] provide simple solutions to create and operate multiple routing
planes, or ”soft” slices [7], to support various QoS and operational
requirements. Each routing plane minimizes a different additive link metric
(e.g., IGP cost, TE metric, QoS) and can exclude some links when needed (i.e.,
constraints in FlexAlgo or exclusions in MTR). In theory, MTR and FlexAlgo
only require a software upgrade of devices to update the control plane, and
they can still rely on standard IP forwarding in the data plane. However, in
practice, this is only the case for MTR, as FlexAlgo is always shipped along
with SRv6 or SR-MPLS. Indeed, FlexAlgo’s implementations rely on SR labels to
identify routing planes and forward traffic consistently on top of them. When
SR is available, TE policies can be configured at ingress nodes to provide
source-based routing on top of existing IGP routing planes (or “algos” when
FlexAlgo is rolled out). SR policies consist of a label stack that is injected
in packet headers to encode the sequence of shortest-path segments that need
to be traversed by packets. In this case, a data plane upgrade, i.e., device
replacement or microcode update whenever possible, is necessary to support SR
forwarding. Indeed, many legacy devices cannot support SR, as mentioned in
[8].

In this paper, we consider deployment scenarios where most devices still
operate with legacy data planes (i.e., without SR support). In this context,
IGP extensions like MTR can be used to create multiple routing planes, but
most routers have not been replaced by SR-capable devices. The most popular
IGP protocol is Open Shortest Path First (OSPF) [9]. In this protocol, routers
periodically broadcast Link State Advertisement (LSA) messages that contain
metrics related to their adjacent links. Upon reception of LSA messages,
routers store them in an LSDB (Link State DataBase) and maintain a Shortest
Path Tree (SPT) towards all the other nodes, which decides the next hop to use
for all possible destinations. This information is stored in the data plane in
a Forwarding Information Base (FIB) and used to forward packets. When only one
link metric is used (e.g., IGP cost), a single routing plane is available. The
set of link metrics is generally designed [10] to satisfy some traffic
scenarios. However, with a single routing plane for all applications, it is,
in general, not possible to satisfy a heterogeneous set of QoS requirements.

To increase routing options, MTR offers the possibility of executing multiple
IGP instances in parallel, each working over a different metric. The SPT
maintained for each instance, or topology, generally aims at satisfying a set
of QoS requirements (e.g., cost, hop count, delay) and administrative
constraints (e.g., exclusion of links or nodes). FlexAlgo [6] can be seen as
an extension of MTR and adds the possibility to build a routing plane over 1)
a given metric, 2) considering some constraints (i.e., exclusions of
nodes/links), and 3) a particular shortest path algorithm. Each algo is
defined by the tuple (metric, constraints, SPT algorithm). Each MTR topology,
or “algo” in FlexAlgo, has its own IP prefix so that traffic can be forwarded
consistently in a hop-by-hop fashion. In both cases, link metrics also need to
be periodically advertised using LSAs so that routers can update their TE
database. To meet a variety of QoS constraints from applications, a large
number of MTR topologies may be required to satisfy all application
requirements. However, multiplying the sets of TE metrics, i.e., MTR
topologies, yields a significant network overhead as LSA messages are sent for
each metric. On top of that, as network conditions fluctuate over time,
topologies may require frequent updates from the management system to follow
the evolution of traffic and network performance. However, link weights can
only be updated at a slow pace, and they cannot follow the evolution of
network quality over time.

To overcome this issue, we propose a new scalable and versatile IGP extension,
called virtual MTR (vMTR), to dynamically maintain a set of virtual routing
topologies. This solution aims at increasing routing flexibility in multi-
topology routing to follow evolving network conditions and support
heterogeneous QoS requirements while reducing the overhead from LSAs and
weight updates. The main idea is to route flows over virtual topologies that
are silent, i.e., they do not exchange LSA messages, and that are continuously
derived from a very limited set of real MTR topologies, i.e., exchanging LSAs.
Each real topologies computes shortest paths for a basic QoS parameter, e.g.
delay or packet loss, while virtual ones are designed to generate paths
meeting multiple QoS constraints.

Searching for a shortest path satisfying one or multiple end-to-end
constraints such as latency or hop count usually involves algorithms relying
on Lagrangian relaxation. These algorithms aim at finding feasible Lagrangian
multipliers to modify link costs so that a shortest path provides a feasible
solution. Well-known heuristic algorithms are LARAC [11], for one constraint,
or GEN-LARAC [12] for multiple constraints. Inspired by this, we propose to
derive virtual topologies from real ones applying a linear combination of
their link metrics with designed multipliers. Compared to vanilla MTR, our
vMTR solution reduces the overhead as each virtual topology is only defined by
a set of multipliers, one per real IGP instance, and virtual topologies do not
exchange LSAs. Furthermore, routing is automatically adjusted when basic QoS
metrics vary and specific link metrics do not need to be updated when traffic
conditions evolve.

In this paper, we first formulate the optimization problem to find optimal
multipliers so as to minimize the number of virtual topologies in vMTR. We
propose a polynomial and exact algorithm that finds the minimum set of
Lagrangian multipliers to meet two QoS constraints for a set of demands. Then,
as a benchmark, we revisit the IGP weight design problem [10] to minimize the
number of topologies in MTR and derive a heuristic algorithm based on local
search to solve the problem. Finally, we compare MTR and vMTR with numerical
results and show that the number of real MTR topologies can be drastically
reduced with vMTR when accepting all traffic demands. We also observe that
vMTR topologies are more robust to QoS changes than MTR ones.

The rest of this paper is structured as follows. Related work is discussed in
Sec. II. The system architecture and the problem formulation are introduced in
Sec. III and IV. Sec. V presents algorithmic solutions. Sec. VI describes our
numerical results and Sec. VII concludes this paper.

##  II Related work

The design of link weights so that traffic requirements can be satisfied is
known as the IGP weight design problem. This problem is a well-studied NP-hard
problem [13].

Single topology routing. Fortz and Thorup [10, 14] have initially proposed a
local search algorithm where the main idea is to change a few link weights at
each time step. They minimize a piece-wise linear cost function related to MLU
(Maximum Link Utilization). Right after, a genetic algorithm has been proposed
by Ericsson et al.[15] and later combined [16] with the original local search.
Bley et al. [13] proposed an improved version of the local search algorithm. A
math-heuristic based on column generation has also been proposed recently
[17].

Other variants related to energy-minimization have been considered [18, 19].
Robust optimization variants have been studied to consider (weighted) link
failures [20] or worst case variations of demands [21] (i.e., minimum /
maximum values for each demand) and hoses [22] at ingress routers to bound
aggregate traffic. Oblivious routing extensions where traffic can be load
balanced over multiple paths have also been considered [23, 24] using the hose
model. Other traffic polytopes have been used to model traffic uncertainty
[25].

Beyond traditional optimization techniques, machine learning [26] have
recently been used to solve the deterministic version of the problem. It is
based on a smoothed formulation of the problem that is suitable for gradient
descent.

Multi-topology routing. In the context of network slicing where end-to-end
latency constraints must be met, a global optimization model along with
results obtained with a genetic algorithm has been presented [7] to design
multiple topologies (iteratively) and assign demands to them. Some papers have
also proposed solution where topologies are designed offline to maximize edge-
disjointness [27, 28] or path diversity [29]. These topologies are then used
to route traffic based on actual QoS conditions. Other variants for MTR have
been studied to consider robustness to traffic variations [30] and failure
scenarios [31].

Our contributions. In this paper, we extend MTR with virtual topologies that
are silent, i.e., they do not exchange LSA messages, and that are continuously
derived from a very limited set of real IGP topologies, i.e., exchanging LSAs.

##  III MTR and vMTR for QoS Constraints

This section introduces MTR, and the vMTR extension we propose, from a system
perspective. As a preliminary step, we also first formulate the QoS routing
problem they address.

Symbol |  Meaning
---|---
$K$ |  Set of demands
$s_{k}$, $d_{k}$ |  Source and destination for demand $k\in K$
$r_{a}^{t}$ |  Resource on arc $a\in A$ for topology $t$
$T_{Q}$ |  Set of basic topologies (each related to a different QoS metric)
$T_{D}$ |  Set of designed topologies
$r_{k}^{t}$ |  Max resource consumption for demand $k$ over topology $t$
$SPT^{t}$ |  Shortest path tree on topology $t$
$\lambda_{t}$ |  Vector of $|T_{Q}|$ multipliers to design vMTR topology $t$
TABLE I: List of key notations.

###  III-A QoS satisfaction problem

Let’s represent a network as a directed graph $G=(V,A)$ where $V$ is the set
of nodes (i.e., routers) and $A$ is the set of arcs (i.e., links). A set of
basic IGP topologies $T_{Q}$ is associated with different QoS metrics such as
delay, jitter or losses. On top of them, additional topologies, real for MTR
or virtual for vMTR, can be designed to route demands that cannot be accepted
on basic topologies. In this case, these topologies defined by the set of
their link weights, also called resources, denoted $r_{a}^{t}$ for arc $a\in
A$ and topology $t\in T_{Q}$, need to be designed along with their mapping to
demands. The set of traffic demands $K$ is given as input. Each demand $k$
goes from source node $s_{k}$ to destination node $d_{k}$, and it has a set of
QoS requirements that are defined by a maximum resource consumption
$r_{k}^{t}$ over basic (additive) QoS metrics from topology $t\in T_{Q}$. Note
that packet loss can be made additive using the log function. The problem is
to find, for every demand $k$, a routing path $p_{k}$ satisfying a set of QoS
constraints, as expressed as follows:

$\displaystyle\sum_{a\in p_{k}}r_{a}^{t}\leq r_{k}^{t}$ $\displaystyle\forall
t\in T_{Q}$ (1)

We consider a practical use case scenario where the total amount of traffic
with QoS constraints is much smaller than best effort traffic. QoS constrained
traffic is prioritized as high priority traffic in the data plane and the sum
of high priority traffic (i.e., QoS constrained traffic) is not enough to
saturate the bottleneck capacity. For this reason, we omit capacity
constraints in the design of topologies and we focus only on QoS constraints
satisfaction.

###  III-B MTR: Multi-Topology Routing

When relying on multi-topology routing, the goal is to design the minimum
number of topologies so that all demands have their QoS requirements
satisfied. As mentioned in Sec. II, this is typically realized by solving a
series of inverse shortest path trees problems. When link weights are designed
and deployed, traffic is routed over a set of $T_{D}$ topologies, each one
corresponding to a new IGP instance. Routers exchange periodically LSA
messages to propagate weights and they maintain a minimum spanning tree
$SPT^{t}$ for every $t\in T_{D}$. Whenever network conditions evolve and the
characteristic of links in terms of QoS changes, a significant update of the
designed topologies, i.e. weights, may be required.

###  III-C vMTR: virtual MTR

In the case of vMTR, additional virtual topologies need to be designed from
the set of existing real topologies. The central idea behind vMTR is rooted in
Lagrangian relaxation, similarly to the LARAC [11] or GEN-LARAC [12]
algorithms that find a path satisfying one or multiple QoS constraints,
respectively. Indeed, virtual topologies are derived from real ones using a
set of Lagrangian multipliers. A set of multipliers, each one associated with
a real topology, is used to calculate and maintain link weights in a virtual
topology. In this way, virtual topologies are automatically updated when QoS
metrics on real instances change. In addition, virtual topologies are
completely silent and do not exchange LSAs.

For instance, a virtual topology can be created from two real topologies on
loss and delay so that (virtual) link weights are defined by $\lambda^{l}$.
Loss + $\lambda^{d}$. Delay with a vector $\lambda$ of multipliers associated
to the topology. Virtual link weights are used by the shortest path tree
algorithm inside routers to find constrained paths with end-to-end loss and
delay constraints. Recall that, in practice, packet loss is made additive
using the log function. In the general case, metrics from real topologies can
be combined with any algebraic operation. However, without loss of generality,
we will consider in the rest of the paper a linear combination. Considering
other operations may help to further increase routing diversity and constraint
satisfaction, but at the cost of some extra complexity. Indeed, Lagrangian
algorithms are actually heuristic to solve the constrained shortest path
problem and may not find feasible solutions when they exist (see Sec. VI).
However, by design, vMTR based on a linear combination is guaranteed to route
all demands for which paths can be found with LARAC.

Figure 1: System overview for vMTR.

To create and manage virtual topologies, a virtual weight designer module sits
in the network controller, as depicted in Fig. 1. It takes as inputs the
network topology, link metrics related to “real” IGP instances (e.g., cost,
capacity, delay), a traffic matrix with a set of demands and their attributes
(e.g., source and, destination, QoS requirements). It produces a set of
$|T_{D}|$ vectors of multipliers, one for every virtual topologies, along with
the mapping of demands to each topology. In practice it sends to all routers
$|T_{D}|$ vectors to create $|T_{D}|$ virtual topologies, each vector being of
size $|T_{Q}|$, i.e., the number of existing / real IGP topologies. Every
router continuously maintains a LSDB with the set of virtual link weights
derived from real metrics. As virtual topologies are based on real ones, no
specific LSAs need to be advertised. These topologies are almost immediately
available to route traffic as routers do not need to wait for LSA exchanges.
When metrics in real topologies are updated, all the virtual weights and their
associated shortest path trees are updated locally by devices. Traditional
Dynamic SPT [32] mechanisms implemented in IGPs can be used to stabilize the
system in case of frequent updates. No routing loops can occur as multipliers
and virtual metrics are positive. Unlike link weights in real topologies that
need to be exchanged with LSAs, multipliers and the resulting virtual weights
do not need to be integer. As for MTR, packets can be routed over a specific
topology using QoS header fields like DSCP or ToS.

In terms of benefits and impact, one advantage of vMRT is the capability to
explain a virtual topology as it is a combination of basic topologies related
to metrics that users understand. This is impossible with MTR when weights are
optimized by an algorithm. Actually, vMTR offers the possibility of manually
or automatically setting multipliers to define a good tradeoff among basic
topologies. An expert or customer can select the best multipliers according to
his requirements. With regards to FlexAlgo, vMTR can be seen as a way to
create new topologies (i.e., new metrics) over which ”algos” are applied. We
could even imagine an extension to the standard where multipliers are part of
the algorithm definition and can be use to tune ”algos” to meet QoS
constraints. However, in this case, SR has to be rolled out.

Our proposition is in this paper two-fold: vMTR provides a way to create
virtual topologies using multipliers and, on top of that, we propose in Sec.
V-A an algorithm to automatically generate the said multipliers.

##  IV Topology design problems

The optimization problem associated with the design of topologies for MTR,
defined by link weights, consists in creating the minimum number of topologies
such that all demands are accepted. This problem is similar in vMTR but
differs in the way weights are computed. Indeed, weights of virtual topologies
are derived from real topologies. This section presents the ILP (Integer
Linear Programming) models for both cases.

In the rest of the paper, we will assume that MTR and vMTR are used to route
demands that cannot be accepted over basic topologies $|T_{Q}|$ (each one
maintaining a shortest path tree over a single QoS metric). Also, to ease the
presentation, we will omit demands for which a constrained path cannot be
found.

$\displaystyle\min\quad$ $\displaystyle\sum_{t\in T_{D}}z_{t}$ (2)
$\displaystyle\sum_{k\in K}y_{k,t}\leq z_{t}$ $\displaystyle\forall t\in
T_{D}$ (3) $\displaystyle\sum_{t\in T_{D}}y_{k,t}\geq 1$ $\displaystyle\forall
k\in K$ (4)
$\displaystyle\sum_{a\in\omega^{+}(v)}x_{a,k,t}-\sum_{a\in\omega^{-}(v)}x_{a,k,t}=\begin{cases}y_{k,t}\text{
if }v=s_{k},\\\ -y_{k,t}\text{ if }v=d_{k},\\\ 0\text{ otherwise;}\end{cases}$
$\displaystyle\forall v\in V,\forall k\in K,\forall t\in T_{D}$ (5)
$\displaystyle x_{a,k,t}\leq u_{a,d_{k},t}$ $\displaystyle\forall a\in
A,\forall k\in K,\forall t\in T_{D}$ (6)
$\displaystyle\sum_{a\in\omega^{+}(v)}u_{a,d,t}\leq 1$ $\displaystyle\forall
v\in V,\forall d\in d_{K},\forall t\in T_{D}$ (7) $\displaystyle
u_{a,d,t}\leq\sum_{k\in K:d=d_{k}}x_{a,d_{k},t}$ $\displaystyle\forall a\in
A,\forall d\in d_{K},\forall t\in T_{D}$ (8) $\displaystyle
w_{(u,v),t}-\pi_{u,d,t}+\pi_{v,d,t}\geq 1-u_{(u,v),d,t}$
$\displaystyle\forall(u,v)\in A,\forall d\in d_{K},\forall t\in T_{D}$ (9)
$\displaystyle w_{(u,v),t}-\pi_{u,d,t}+\pi_{v,d,t}\leq M(1-u_{(u,v),d,t})$
$\displaystyle\forall(u,v)\in A,\forall d\in d_{K},\forall t\in T_{D}$ (10)
$\displaystyle\sum_{a\in A}\sum_{t\in T_{D}}r_{a}^{t}x_{a,k,t}\leq r_{k}^{t}$
$\displaystyle\forall k\in K,\forall t\in T_{Q}$ (11) Figure 2: Integer linear
problem for the design of topologies in MTR.

###  IV-A MTR

Let us introduce the variables to solve the topology design problem for MTR:

  * •

$x_{a,k,t}\in\\{0,1\\}$: equals 1 if the demand $k$ is routed along a path
using link $a$ in the topology $t$, and 0 otherwise.

  * •

$\pi_{u,v,t}\in\mathbb{R}^{+}$: is the potential of the node $u$, i.e., the
distance between node $u$ and $v$ on the topology $t$.

  * •

$u_{a,d,t}\in\\{0,1\\}$: equals 1 if the link $a$ belongs to a shortest path
towards destination $d$ on topology $t$, and 0 otherwise.

  * •

$w_{a,t}\in\mathbb{R}^{+}$: the weight assigned to link $a\in A$ on the
topology $t$.

  * •

$y_{k,t}\in\\{0,1\\}$: equals 1 if the demand $k$ is assigned to the topology
$t$ (thus feasible for this topology), and 0 otherwise.

  * •

$z_{t}\in\\{0,1\\}$: equals 1 if the topology $t$ is deployed.

Let us consider $d_{K}$ the set of nodes that are destinations of at least one
demand. We denote by $\omega^{-}(v)$ (resp. $\omega^{+}(v)$) the set of
ingoing (resp. outgoing) links of node $v$.

As we do not know the needed number of topologies, we can consider an upper
bound to create the right number of variables. This upper bound can be set by
the number of designed topologies returned by the local search algorithm
presented in the next section.

The ILP model is given by Fig. 2. This model is an extension of the well-known
model from the literature (e.g., [30, 7]) to jointly create multiple
topologies. The main goal is to minimize the number of designed topologies.
Constraints (3) link the activation of topologies and demands routed on top of
them. Inequalities (4) ensure that each demand is assigned to one designed
topology. Equalities (5) are the well-known flow conservation constraints.
Constraints (6) and (8) link routing variables and shortest path variables.
Inequalities (7) guarantee that the induced shortest path tree is a tree,
i.e., there is at most one outgoing link at each node. Constraints (9) and
(10) ensure that the shortest path tree follows the weight $w$. Inequalities
(12) are resource constraints.

###  IV-B vMTR

For vMTR we need to add the following variables:

  * •

$\lambda_{t,\bar{t}}\in\mathbb{R}$: the multiplier associated with the basic
topology $t\in T_{Q}$ and the virtual topology $\bar{t}\in T_{D}$.

and the following constraints:

$\displaystyle w_{a,\bar{t}}=\sum_{t\in T_{Q}}\lambda_{t,\bar{t}}w_{a,t}$
$\displaystyle\forall a\in A,\forall\bar{t}\in T_{D}$ (12)

Remark that in the case of vMTR, we complete the designed virtual topologies
with real MTR topologies when all demands cannot be accepted with vMTR. To
manage the two kinds of designed topologies, we penalize in the objective
function the activation of MTR topologies.

##  V Algorithmic solutions

This section presents algorithms for MTR and vMTR. While multiple algorithms
exist for MTR (see Sec. II), we build upon one of them [13] to focus on QoS
constraint satisfaction and minimize the number of designed topologies.

###  V-A MTR Algorithm

Algorithm 1 Greedy algorithm for MTR

1:function greedy($G=(V,A),d,l,K$)

2: $\mathcal{T}\leftarrow\emptyset$ $\triangleright$ Set of weight-demands
pairs

3: while $K\neq\emptyset$ do

4: $w\leftarrow\\{\mathcal{U}(0,65535)|a\in A\\}$$\triangleright$ random
values

5: $w\leftarrow\textsc{localSearch}(G,w,d,l,K)$

6: $K^{\prime}\leftarrow\textsc{getAcceptedDemands}(G,w,d,l,K)$

7: if $K^{\prime}=\emptyset$ then

8: $p\leftarrow\textsc{TAMCRA}(G,d,l,K)$

9: $w\leftarrow\\{\mathcal{U}(0,65535)|a\in A\\}$ $\triangleright$ random
values

10: for all $a\in p$ do

11: $w_{a}\leftarrow 1$

12: $K^{\prime}\leftarrow\textsc{getAcceptedDemands}(G,w,d,l,K)$

13: $\mathcal{T}\in\mathcal{T}\cup{(w,K^{\prime})}$

14: $K\leftarrow K\setminus K^{\prime}$

15:function localSearch($G=(V,A),w,d,l,K$)

16: $\tilde{K}\leftarrow\textsc{getAcceptedDemands}(G,w,d,l,K)$

17: $(w^{\ast},K^{\ast})\leftarrow(w,\textsc{getAcceptedDemands}(G,w,d,l,K))$

18: for $i\leftarrow 0$ to maxIte do

19: $\mathcal{W}\leftarrow\textsc{generateNeighborhood}(w,G,d,l)$

20:
$w\leftarrow\arg\max_{w\in\mathcal{W}}\textsc{getAcceptedDemands}(G,w,d,l)$

21: $\tilde{K}\leftarrow$getAcceptedDemands($G,w,d,l,K$)

22: if $|\tilde{K}|>|K^{\ast}|$ then

23: $(w^{\ast},K^{\ast})\leftarrow(w,\tilde{K})$

24: return $w^{\ast}$

25:function generateNeighborhood($w,G=(V,A)$)

26: $(\delta^{+},\delta^{-})\leftarrow$computeDeltaWeights(G, w)

27: $\mathcal{W}\leftarrow\emptyset$

28: for all $a\in A$ do

29: $w^{-}\leftarrow w$

30: $w_{a}^{-}\leftarrow w_{a}+\delta^{-}(a$)

31: $w^{+}\leftarrow w$

32: $w_{a}^{+}\leftarrow w_{a}+\delta^{+}(a)$

33: $\mathcal{W}\leftarrow\mathcal{W}\cup\\{w^{+},w^{-}\\}$

34: return $\mathcal{W}$

35:function computeDeltaWeights($G=(V,A),w$)

36: $\delta^{+}\leftarrow\\{\infty|a\in A\\}$

37: $\delta^{-}\leftarrow\\{\infty|a\in A\\}$

38: for all $v\in V$ do

39: $T,d\leftarrow\textsc{dijkstraShortestPathTree}(G,v,w)$ $\triangleright$
Get shortest path tree with distances

40: for all $(s_{a},t_{a})\in A\setminus T$ do

41: $v_{A}\leftarrow\textsc{LowestCommonAncestor}(T,s_{a},t_{a})$

42: for all $a^{\prime}\in p_{T(v_{A},s_{a})}$ do

43:
$\delta^{-}(a^{\prime})\leftarrow\min(\delta^{-}(a^{\prime}),d(s_{a})+w(a)-d(t_{a}))$

44: for all $a^{\prime}\in p_{T(v_{A},t_{a})}$ do

45:
$\delta^{+}(a^{\prime})\leftarrow\min(\delta^{+}(a^{\prime}),d(s_{a})+w(a)-d(t_{a}))$

We developed a greedy algorithm to design the sets of weights, i.e.,
topologies, as shown in Algo. 1. It iteratively generates random weight
vectors, fine-tuned to maximize the number of accepted demands using a local
search algorithm. MTR topologies are created one by one, each time trying to
accept a maximum number of demands. If the local search fails to generate a
weight vector accepting any demands, it falls back to computing a valid
resource-constrained shortest path for one of the remaining demands (using the
TAMCRA algorithm [33]), then uses that path for building a valid weight vector
from a new random one.

The greedy algorithm uses a run-of-the-mill local search that encodes a
solution by its weight vector, and the neighborhood of a solution is the set
of weight vectors that differ from the solution by a single element. The
fitness function used to evaluate a solution computes the shortest path tree
according to the weight vector and counts the number of demands with a path
respecting their required QoS.

Bley et al. [13] already proposed a local search algorithm for the IGP weight
design problem, but we consider different assumptions. Firstly, while they
consider an extended neighborhood (compared to ours) that takes into account
load balancing, we only focus on the weight vector variation. Secondly, while
they consider weights in the range $[0,20]$, we consider the whole range of
possible OSPF weights ($[0,2^{16}]$) The last point is particularly important
since the definition of their neighborhood considers atomic changes in the
weight ($\pm 1$), which is quite effective considering the range of the
weights. However, since we allow a much larger range of weights, we need to
define a neighborhood that guarantees changes in the shortest path trees
induced by the weights.

We thus generate meaningful neighbors (see line 25) computing the minimum
increase or decrease in weight needed to change the shortest path trees. We
call these values the up and the down delta-weight of a link, respectively. We
denote by $\pi(u,v)$ the distance of the shortest path from $u$ to $v$ with
the current link weights.

First, we define the down delta-weight $\delta^{-}_{v}(a)$ of a link
$a=(s_{a},d_{a})$ as the weight decrease needed to insert the link $a$ in the
current shortest path tree $T_{v}$, rooted at ${v}$, in the place of the link
$a^{\prime}=(s_{a}^{\prime},d_{a})$. Formally, we are looking for a value
$\delta^{-}_{v}(a)$ such that

$\pi(v,s_{a})+w_{(a)}-\delta^{-}_{v}(a)<\pi(v,t_{a}),$
$\delta^{-}_{v}(a)=\pi(v,s_{a})+w_{(a)}-\pi(v,t_{a})$

Second, we define the up delta-weight $\delta^{+}_{v}(a)$ of a link
$a=(s_{a},d_{a})$ as the weight increase needed to remove the link from the
current shortest path tree $SPT_{v}$, rooted at $v$. Formally, we are looking
for a value $\delta^{+}_{v}(a)$ such that

$\pi(v,s_{a})+w_{(a)}+\delta^{+}_{v}(a)>\min_{a^{\prime}\in\omega^{-}(d_{a})\setminus\\{a\\}}(\pi(v,s_{a^{\prime}})+w_{(a^{\prime})}).$

Thus, the up delta weight of a link $a$ is given by

$\displaystyle\delta^{+}_{v}(a)$
$\displaystyle=\min_{a^{\prime}\in\omega^{-}(d_{a})\setminus\\{a\\}}(\pi(v,s_{a^{\prime}})+w_{(a^{\prime})})-\pi(v,s_{a})-w_{(a)}$
(13)
$\displaystyle=\min_{a^{\prime}\in\omega^{-}(d_{a})\setminus\\{a\\}}\delta^{-}_{v}(a^{\prime})$
(14)

Once we computed the delta-weight of each link, we can propagate them (see
line 40) up the tree. Indeed, any weight change occurring in the higher depth
of the tree will impact the distance of the descendants. Given a link
$a_{1}=(s_{1},d_{a_{1}})$ and its swapped link $a_{2}=(s_{2},d_{a_{1}})$, we
can propagate their up and down delta-weight up to their lowest common
ancestor in the tree. Formally, if $u$ is the lowest common ancestor of
$s_{a_{1}}$ and $s_{a_{2}}$ in $SPT$, and $p_{T}(v_{1},v_{2})$ is the path
from $v_{1}$ to $v_{2}$ in $SPT$, we have the following properties:

$\forall a\in p_{T}(u,s_{1}),\delta^{+}(a)=\delta^{+}(a_{1}),$ $\forall a\in
p_{T}(u,s_{2}),\delta^{-}(a)=\delta^{-}(a_{2})$

Finally, we need to consider all possible shortest path trees in the graph. We
can combine their up and down delta-weight by taking the minimum over all
trees.

###  V-B vMTR Algorithm

The weight design problem for vMTR consists in maximizing the number of
accepted demands (first objective) while minimizing the number of virtual
topologies (second objective). To create a virtual topology with vMTR over 2
basic topologies (e.g., delay and loss), we consider a single multiplier
$\lambda$ that provides the following virtual weight
$r_{a}^{t}=r_{a}^{1}+\lambda r_{a}^{2}$ for each link $a\in A$. A multiplier
$\lambda$ is said to be feasible for demand $k$ if the shortest path $p$ over
link metrics $r_{a}^{t}$ satisfies the maximum constraint on the end-to-end
consumption for both resources, i.e., on $r_{k}^{1}$ and $r_{k}^{2}$.

In the following, we first present some properties for this particular case
with 2 resources and then, we propose an exact algorithm that solves the
problem in polynomial time.

Lagrangian properties. As one multiplier can cover several demands, a natural
question is: does there exist one or several intervals of feasible multipliers
for a given demand? In this section, we provide an answer to this question
that allows solving the vMTR problem in polynomial time afterward.

Let us introduce the Constraint Shortest Path problem (CSP) where the goal is
to find a path minimizing a resource under an additive end-to-end resource
constraint. In [34], the authors propose and analyze the LARAC algorithm, a
polynomial time heuristic solution for the CSP problem based on Lagrangian
relaxation. This algorithm computes a multiplier $\lambda^{*}$ and the
solution it returns is the shortest path $p_{\lambda^{*}}$ in the graph where
the cost of each link $a\in A$ is given by $r_{a}^{1}+\lambda^{*}r_{a}^{2}$.

The most interesting claim proved in [34] for us is given below.

######  Claim 1.

[34] If $\lambda<\lambda^{*}$, then $\sum_{a\in p_{\lambda}}r^{2}_{a}\geq
r^{2}_{k}$ and if $\lambda>\lambda^{*}$, then $\sum_{a\in
p_{\lambda}}r^{2}_{a}\leq r^{2}$

This claim is the core of our vMTR algorithm.

We denote by $\lambda^{1*}_{k}$ (resp. $\lambda^{2*}_{k}$) the multiplier
provided by LARAC algorithm when resource 2 (resp. resource 1) is minimized
and resource 1 (resp. resource 2) is constrained. As proved in [34] this
multiplier is optimal.

The following proposition ensures that, for a given demand, only one range of
multiplier is feasible for the two resources.

######  Proposition 1.

For a given demand $k\in K$, a multiplier $\lambda_{k}$ is feasible for $k$,
if

$\lambda^{2*}_{k}\leq\lambda_{k}\leq\frac{1}{\lambda^{1*}_{k}}.$

###### Proof.

Let us consider the shortest path $p$ where link costs
$r_{a}^{2}+\lambda^{1*}_{k}r_{a}^{1},\forall a\in A$ are considered. The path
$p$ is also the shortest path in the graph where link costs are
$r_{a}^{1}+\frac{1}{\lambda^{1*}_{k}}r_{a}^{2},\forall a\in A$. Indeed, for
each $a\in A$, dividing $r_{a}^{2}+\lambda^{1*}_{k}r_{a}^{1}$ by
$\lambda^{1*}_{k}$ induced $r_{a}^{1}+\frac{1}{\lambda^{1*}_{k}}r_{a}^{2}$.
Thus, thanks to the claim 1, all $\lambda^{2}_{k}\geq\lambda^{2*}_{k}$ induced
a shortest path where resource 2 is respected and all
$\frac{1}{\lambda^{1}_{k}}\leq\frac{1}{\lambda^{1*}_{k}}$ induced a shortest
path where resource 1 is respected. Thus all $\lambda_{k}$ such that
$\lambda^{2*}_{k}\leq\lambda_{k}\leq\frac{1}{\lambda^{1*}_{k}}$ are feasible
for the two resources. ∎

For a given multiplier $\lambda$, we may have several paths that are the
shortest, but one of them can be unfeasible for resource 1 or 2. To address
this issue, we modified the Dijkstra algorithm to find all shortest paths
between two nodes in a graph and we slightly changed the $\lambda$ value by
adding an $\epsilon$ until all shortest paths are valid or only one path is
the shortest.

Figure 3: Example of non shortest paths uniqueness with $\lambda=13$.

In the following, for a given demand $k$, we denote
$\lambda_{k}^{\min}=\lambda^{2*}_{k}$ and
$\lambda_{k}^{\max}=\frac{1}{\lambda^{1*}_{k}}$. Thanks to the convergence of
LARAC algorithm [34], $\lambda_{k}^{\min}$ is minimum for each demand $k\in K$
and $\lambda_{k}^{\max}$ is maximum for each demand $k\in K$.

From Proposition 1, we can derive an efficient polynomial algorithm to solve
the weight design problem for vMTR. As a pre-processing step, we compute the
multipliers $\lambda_{k}^{\min}$ and $\lambda_{k}^{\max}$ for each demand
$k\in K$. Each demand where $\lambda_{k}^{\min}>\lambda_{k}^{\max}$ is
discarded as it cannot be routed with vMTR. At the end of the procedure, the
pool of discarded demand is processed with the MTR algorithm (Sec. V-A). The
overall goal is to minimize the number of MTR topologies, as they generate
LSAs and do not adapt to evolving network conditions.

Exact algorithm. The goal of this algorithm is to find a set of multipliers
$\Lambda$ such that all demands are covered. A multiplier $\lambda$ covers a
demand $k$ if $\lambda$ is feasible for the demand (see Proposition 1). At
each step, we find a new multiplier to add in $\Lambda$ such that some
uncovered demands become covered with the new multiplier. More formally. Let
$\bar{K}$ be the set of demands covered at the current step, and respectively,
$K\setminus\bar{K}$ the set of uncovered demands.

Algo. 2 describes the whole algorithm and an example is given in Fig. 4.
First, the algorithm orders the $\lambda_{k}^{\max}$ in increasing order.
Then, at each step, it selects a $\lambda_{k^{*}}^{\max}$ to add in $\Lambda$
such that all uncovered commodities have an associated $\lambda_{k}^{\max}$
greater or equal to the selected $\lambda_{k^{*}}^{\max}$. It continues until
all demands are covered.

Figure 4: Example with $6$ demands. Intervals between $\lambda_{k}^{\min}$ and
$\lambda_{k}^{\max}$ are represented for each demand $k$. The first multiplier
$\lambda_{1}$ covers the green demands, the second multiplier $\lambda_{2}$
covers the blue demands and $\lambda_{3}$ covers the orange demand. Algorithm
2 Weight design algorithm for vMTR

1:function ($G=(V,A),K$)

2: $\bar{K}=\emptyset$, $\Lambda=\emptyset$

3: $\sigma$ ordering of demands $k\in K$ in increasing order of
$\lambda_{k}^{\max}$

4: while $K\setminus\bar{K}\neq\emptyset$ do

5: $\lambda^{*}$ the smallest $\lambda_{k}^{\max}$ such that $k\in
K\setminus\bar{K}$

6: add in $\bar{K}$ all demands covered by $\lambda^{*}$ not in $\bar{K}$

7: add $\lambda^{*}$ in $\Lambda$

8: return $\Lambda$

Complexity and generalization. To summarize the approach above, the problem is
solved in two steps. First, we define the space of valid multipliers for each
demand. Second, we search the minimum number of multiplier vectors such that
all demands are covered.

The dimension of the space of valid multipliers is equal to $|T_{Q}|-1$. For
the case with $2$ basic topologies, we have only 1 dimension, which can be
represented by an interval graph and the problem is polynomial as presented
above. If we consider $3$ basic topologies, we have 2 dimensions that can be
represented by a plane. In this case, we can consider that the area of valid
multipliers, vector of size 2, represents a rectangle, which is a particular
case of our problem. This kind of graph is called boxicity 2-graph. For the
interval graph or a boxicity 2-graph, the selection of the minimum number of
multiplier vectors consists in solving a clique cover problem in the graph
$H=(K,I)$ where $K$ is the set of demands and $I$ is the set of links. Two
demands are linked if and only if they share a common valid multiplier vector.
This problem is NP-hard for a boxicity 2-graph which is a particular case of
the problem when 3 basic topologies are considered and polynomial for an
interval graph [35].

##  VI Performance evaluation

We now compare MTR and vMTR over a set of instances generated from SNDlib
[36].

Description of instances. These instances contain the geographical coordinates
for each node and the link capacities between them. We selected 15 large
topologies. For each network, we have considered two resources, the loss and
the delay. For each link, the packet loss is set to be inversely proportional
to the link capacity and the delay is considered to be proportional to the
distance. To avoid trivial instances, we generated all demands between each
pair of routers. We set the maximum end-to-end loss (resp. delay) to one of
the shortest path on the loss (resp. delay) minus $\epsilon$. This ensures
that the demands cannot be routed on basic topologies, i.e., loss or delay.
For each demand, we solved optimally, i.e., not using LARAC, the constraint
shortest path problem to ensure that a solution exists. If no solution exists
for a demand then we did not consider it. All instances are available at [37].

[width=]figures/nb_topologies

(a) Number of topologies.

[width=]figures/average_demands_per_topo

(b) Average demands per topology.

[width=]figures/execution_time

(c) Execution times.

Figure 5: Comparison of MTR and vMTR using violin plots.

Numerical results. Fig. 5 presents violin plots to summarize results over all
15 instances. We can see, in Fig. 5(a), that vMTR reduces the number of real
topologies compared to MTR, from an average of 8.41 to 5.41. It also reduces
the maximum number of needed topologies from 20 to 12, and it also avoid the
need for any real topology in some cases. vMTR requires an average of 4.08
virtual topologies. The decrease in the number of real topologies can be
easily explained by looking at the number of average demands per topologies in
Fig. 5(b). vMTR virtual topologies can accept more demands than MTR real ones:
in average, a vMTR topology holds 32.7 demands while, in average, MTR
topologies hold 18.7 demands. In some cases, we can even have virtual
topologies with 131.4 demands. However, real topologies in vMTR hold less
demands because the remaining demands not covered by a virtual topology are
harder to cover with the greedy MTR algorithm. Recall that as our
implementation of vMTR with a linear combination of real topologies is based
on the same principle as LARAC, which is a heuristic solution to the CSP
problem, some demands cannot be satisfied with virtual topologies and real
ones need to be added to ensure 100% traffic acceptance.

Fig. 5(c) shows the execution times for both solutions. For vMTR, we present
the time required to compute lambda multipliers for virtual topologies as well
as the time needed by the greedy algorithm to provision the remaining demands
over real topologies. vMTR, once again, shows an improvement over MTR with an
average time of 86.2s over 141.3s, and it also decreases the maximum time from
740.5s to 411.1s. The computational time to compute only the virtual
topologies is really small. It is due to the polynomiality of the algorithm.

[width=]figures/qos-robustness

Figure 6: Violin plots of the QoS robustness which is, for a demand and a
metric, the ratio between the end-to-end QoS of the provided paths and the
requested QoS. A value of 1 for a given metric means that the path of a demand
is at its maximum resource consumption.

Finally, we study the robustness of the two approaches against changes of QoS
metrics. To this end, we define the QoS robustness of a demand, for a given
QoS metric, as the ratio of its provisioned path QoS over its end-to-end QoS
requirement. A QoS robustness of 1 means that the provided path QoS equals the
required QoS, and any increase of this QoS metric in the path’s link would
invalidate the path. While we can see in Fig. 6 that both solutions provision
demands on strict paths (ratio of 1), on average, vMTR provides paths that are
more robust to QoS changes than MTR (0.75 and 0.79 for vMTR vs 0.85 and 0.88
for MTR, for delay and losses, respectively). Remark that the robustness is
slightly better from the delay point of view if we compare with the Loss. This
is due to the distribution of metric values in our instances. For the two
solutions, we have a better maximum robustness for the delay than for the
loss. We then have more flexibility to improve the robustness of the delay.

##  VII Conclusion

We introduced a new IGP extension, called virtual MTR (vMTR), to dynamically
maintain a set of virtual routing topologies. This solution aims at increasing
routing flexibility in multi-topology routing to follow evolving network
conditions and support heterogeneous QoS requirements while reducing the
overhead from LSAs and weight updates. We presented a polynomial and exact
algorithm for vMTR when virtual topologies are derived from two real
topologies and a local search algorithm for MTR as benchmark. We showed that
vMTR can help reducing the number of real topologies and is more robust to QoS
changes.

Future work along these lines include several directions. First, we may
propose an efficient algorithm for the design of MTR topologies jointly with
the creation of virtual ones in vMTR. We can also design lambda multipliers
considering traffic scenarios or uncertainty sets on traffic. Also of
interest, we may consider non additive and non linear ways to combine real
topologies. Finally, we may develop efficient algorithms for the hard case
where more than two real topologies are used to derive virtual ones.

## References

  * [1] N. Wang, K. H. Ho, G. Pavlou, and M. Howarth, “An overview of routing optimization for internet traffic engineering,” _IEEE Communications Surveys Tutorials_ , vol. 10, no. 1, pp. 36–56, 2008.
  * [2] P. L. Ventre, S. Salsano, M. Polverini, A. Cianfrani, A. Abdelsalam, C. Filsfils, P. Camarillo, and F. Clad, “Segment routing: a comprehensive survey of research activities, standardization efforts, and implementation results,” _IEEE Communications Surveys & Tutorials_, vol. 23, no. 1, pp. 182–221, 2020.
  * [3] D. Wu and L. Cui, “A comprehensive survey on segment routing traffic engineering,” _Digital Communications and Networks_ , 2022.
  * [4] P. Psenak, S. Mirtorabi, A. Roy, L. Nguyen, and P. Pillay-Esnault, “Multi-topology (MT) routing in OSPF. IETF RRC 4915.” Tech. Rep., 2007.
  * [5] T. Przygienda, N. Shen, and N. Sheth, “M-ISIS: multi topology (MT) routing in Intermediate System to Intermediate Systems (IS-ISs). IETF RFC 5120.” Tech. Rep., 2008.
  * [6] P. Psenak, S. Hegde, C. Filsfils, K. Talaulikar, and A. Gulko, “IGP Flexible Algorithm,” RFC 9350, Feb. 2023.
  * [7] N. Huin, S. Martin, J. m. Leguay, and S. Cai, “Network slicing with multi-topology routing,” in _2021 IFIP Networking Conference (IFIP Networking)_. IEEE, 2021, pp. 1–2.
  * [8] Y. Tian, Z. Wang, X. Yin, X. Shi, Y. Guo, H. Geng, and J. Yang, “Traffic engineering in partially deployed segment routing over ipv6 network with deep reinforcement learning,” _IEEE/ACM Transactions on Networking_ , vol. 28, no. 4, pp. 1573–1586, 2020.
  * [9] J. T. Moy, _OSPF: anatomy of an Internet routing protocol_. Addison-Wesley Professional, 1998.
  * [10] B. Fortz and M. Thorup, “Internet traffic engineering by optimizing OSPF weights,” in _Proc. IEEE INFOCOM_ , vol. 2, 2000, pp. 519–528.
  * [11] A. Juttner, B. Szviatovski, I. Mécs, and Z. Rajkó, “Lagrange relaxation based method for the qos routing problem,” in _Proc. IEEE INFOCOM_ , vol. 2, 2001, pp. 859–868.
  * [12] Y. Xiao, K. Thulasiraman, and G. Xue, “Gen-larac: A generalized approach to the constrained shortest path problem under multiple additive constraints,” in _International Symposium on Algorithms and Computation_. Springer, 2005, pp. 92–105.
  * [13] A. Bley, “Routing and capacity optimization for IP networks,” in _Operations Research Proceedings 2007_. Springer, 2008, pp. 9–16.
  * [14] B. Fortz and M. Thorup, “Optimizing OSPF/IS-IS weights in a changing world,” _IEEE journal on selected areas in communications_ , vol. 20, no. 4, pp. 756–767, 2002.
  * [15] M. Ericsson, M. G. C. Resende, and P. M. Pardalos, “A genetic algorithm for the weight setting problem in OSPF routing ,” _Journal of combinatorial optimization_ , vol. 6, no. 3, pp. 299–333, 2002.
  * [16] L. S. Buriol, M. G. Resende, C. C. Ribeiro, and M. Thorup, “A hybrid genetic algorithm for the weight setting problem in OSPF/IS-IS routing,” _Networks: An International Journal_ , vol. 46, no. 1, pp. 36–56, 2005.
  * [17] A. Benhamiche, M. Chopin, and S. Martin, “Unsplittable shortest path routing: Extended model and matheuristic,” in _Proc. CoDIT_ , 2023.
  * [18] A. Capone, C. Cascone, L. G. Gianoli, and B. Sanso, _OSPF optimization via dynamic network management for green IP networks_. IEEE, 2013.
  * [19] E. Amaldi, A. Capone, and L. G. Gianoli, “Energy-aware IP traffic engineering with shortest path routing,” _Computer Networks_ , vol. 57, no. 6, pp. 1503–1517, 2013.
  * [20] B. Fortz, M. Thorup _et al._ , “Robust optimization of OSPF/IS-IS weights,” in _Proc. INoC_ , vol. 20, 2003, pp. 225–230.
  * [21] J. Vallet and O. Brun, “Online OSPF weights optimization in IP networks,” _Computer Networks_ , vol. 60, pp. 1–12, 2014.
  * [22] R. S. Ranaweera, I. M. Kamrul, and E. Oki, “Preventive start-time optimization of OSPF link weights against link failure for hose model,” in _Asia-Pacific Conference on Communications (APCC)_. IEEE, 2012, pp. 698–702.
  * [23] A. Altin, B. Fortz, and H. Ümit, “Oblivious OSPF routing with weight optimization under polyhedral demand uncertainty,” _Networks_ , vol. 60, no. 2, pp. 132–139, 2012.
  * [24] N. G. Duffield, P. Goyal, A. Greenberg, P. Mishra, K. K. Ramakrishnan, and J. E. van der Merive, “A flexible model for resource management in virtual private networks,” _SIGCOMM Comput. Commun. Rev._ , vol. 29, no. 4, p. 95–108, Aug. 1999.
  * [25] J. Moulierac and T. K. Phan, “Optimizing IGP link weights for energy-efficiency in multi-period traffic matrices,” _Computer Communications_ , vol. 61, pp. 79–89, 2015.
  * [26] M. Kodialam and T. Lakshman, “Network link weight setting: A machine learning based approach,” in _Proc. IEEE INFOCOM_ , 2022, pp. 2048–2057.
  * [27] N. Wang, K.-H. Ho, and G. Pavlou, “Adaptive multi-topology IGP based traffic engineering with near-optimal network performance,” in _IFIP NETWORKING_. Springer, 2008, pp. 654–666.
  * [28] N. Wang, K. H. Ho, and G. Pavlou, “Ample: an adaptive traffic engineering system based on virtual routing topologies,” _IEEE Communications Magazine_ , vol. 50, no. 3, pp. 185–191, 2012.
  * [29] E. Mirzamany, A. Lasebae, and O. Gemikonakli, “An efficient traffic engineering based on multi-topology routing for future internet,” _Computer Networks_ , vol. 70, pp. 170–178, 2014.
  * [30] X. Wang, S. Wang, and L. Li, “Robust traffic engineering using multi-topology routing,” in _GLOBECOM 2009-2009 IEEE Global Telecommunications Conference_. IEEE, 2009, pp. 1–6.
  * [31] M. Menth and R. Martin, “Network resilience through multi-topology routing,” in _DRCN_ , 2005, pp. 271–277.
  * [32] D. Garg _et al._ , “Dynamizing dijkstra: A solution to dynamic shortest path problem through retroactive priority queue,” _Journal of King Saud University-Computer and Information Sciences_ , vol. 33, no. 3, pp. 364–373, 2021\.
  * [33] H. De Neve and P. Van Mieghem, “Tamcra: a tunable accuracy multiple constraints routing algorithm,” _Computer Communications_ , vol. 23, no. 7, pp. 667–679, 2000.
  * [34] X. Ying, T. Krishnaiyan, X. Guoliang, A. Jüttner, and S. Arumugam, “The constrained shortest path problem,” _AKCE International Journal of Graphs and Combinatorics_ , vol. 2, no. 2, pp. 63–86, 2005.
  * [35] H. Imai and T. Asano, “Finding the connected components and a maximum clique of an intersection graph of rectangles in the plane,” _J. Algorithms_ , vol. 4, pp. 310–323, 1983.
  * [36] S. Orlowski, R. Wessäly, M. Pióro, and A. Tomaszewski, “SNDlib 1.0—Survivable network design library,” _Networks: An International Journal_ , vol. 55, no. 3, pp. 276–286, 2010.
  * [37] N. Huin, S. Martin, and J. Leguay, “Loss/delay instances for vIGP,” Jan. 2024\. [Online]. Available: https://doi.org/10.5281/zenodo.10469172