Joao C. Pinto Barros

# Towards the Application of Skewed Detailed Balance in Lattice Gauge Theories

Marina Krstic Marinkovic

###### Abstract

State-of-the-art algorithms in lattice gauge theory typically rely heavily on
detailed balance, which is an instrumental tool to prove the correct
convergence of the Markov Chain Monte Carlo Algorithm. In this work, we
investigate an alternative to detailed balance, skewed detailed balance, and
the prospects of alleviating the topological freezing problem by studying the
one-dimensional $O\left(2\right)$ model.

##  1 Introduction

In a Markov Chain Monte Carlo (MCMC) algorithm, ergodicity and detailed
balance guarantee the convergence to the desired probability distribution.
While ergodicity ensures its existence, detailed balance guarantees that the
distribution is the target one. Thus, convergence does not necessitate
detailed balance. MCMC simulations of Quantum Chromodynamics (QCD), and other
theories, suffer from topological freezing (see e.g. [1, 2, 3, 4, 5]), which
is an instance of critical slowing down. This constitutes an effective
ergodicity breaking. While, typically, the algorithm remains ergodic, the time
necessary to sample all relevant configurations grows increasingly large.
Ergodicity is then effectively broken and, with it, the guarantee that the
obtained results can be trusted.

Applying detailed balance breaking algorithms has been proven advantageous in
various situations111Perhaps the most well-known example of detailed balance
breaking is sequential updating in classical spin systems. When applying the
Metropolis algorithm, for example, one can construct an algorithm that chooses
a random spin to be updated at every step or, alternatively, an algorithm that
transverses the lattice sequentially. The latter breaks detailed balance but,
when done correctly, still converges to the correct distribution and is more
efficient. ([6, 7, 8, 9, 10, 11] to name a few). It can be argued intuitively
that it is expected that algorithms breaking detailed balance can outperform
those that do not. Detailed balance corresponds to an equilibrium condition
where all elementary processes are in equilibrium with their inverse. In MCMC
algorithms this implies that it is possible that for large Monte Carlo times,
the algorithm keeps coming to the same configuration, or remains in a certain
"neighborhood", not exploring fast enough the full configuration space. In
turn, if one allows for detailed balance-breaking processes, one may forbid
the updates to flow in a preferred direction in the configuration space,
exploring it faster. This idea is further explored below.

While the ultimate goal is to address this problem in QCD, $O\left(N+1\right)$
models in $N$ dimensions are an interesting setting to explore. While being
much simpler, they also have non-trivial topological sectors. On the other
side, MCMC algorithms that sample effectively the different topological
sectors are well known, as the Wolff cluster algorithm [12, 13]. Here we will
be focusing on a type of topological action [14] for the 1-d $O\left(2\right)$
model.

##  2 A 1-d $O\left(2\right)$opological Action

The 1-d $O\left(2\right)$ model is a simple model exhibiting different
topological sectors. Throughout the text, we shall always assume periodic
boundary conditions. Under those circumstances and in the continuum, the
action can be regarded as a function that takes functions that map the circle
on itself: $\theta:O\left(2\right)\rightarrow O\left(2\right)$. Finite action
configurations satisfy $\theta\left(x+2\pi\right)=\theta\left(x\right)+2\pi
Q$, for some integer $Q$, where $Q$ represents the topological charge and
corresponds to the winding number of the $O\left(2\right)$ vector for a given
configuration. As it is not possible to continuously deform
$\theta\left(x\right)$ with a given winding $Q$ to a new configuration
$\theta^{\prime}\left(x\right)$ with a different winding $Q^{\prime}\neq Q$,
these configurations can be said to be separated by infinite action barriers.

On the lattice, the standard action takes the form
$S_{\mathrm{std}}=-1/g^{2}\sum_{n=1}^{N}\cos\left(\theta_{n+1}-\theta_{n}\right)$.
While the application of local algorithms, like Metropolis-Hastings or heat-
bath, works in a very similar way as long as the action remains local, some
technical difficulties arise when applying skewed detailed balance to lattice
field theories of continuous variables. As a testing ground, we use a
topological lattice action [14]

$S=\sum_{n=1}^{N}s_{\delta}\left(\theta_{n+1}-\theta_{n}\right),\quad\quad
s_{\delta}\left(\Delta\right)=\left\\{\begin{array}[]{ll}0&\cos\left(\Delta\right)<\delta\\\
s_{0}&\mathrm{otherwise}\end{array}\right..$ (1)

It strips down the 1-d $O\left(2\right)$ model to the basic features we want
to explore. The parameter $\delta$ plays the role of the coupling ($g$ in the
standard action). In turn, $s_{0}$ is the height of the barriers that separate
different topological sectors. For a system of size $L$, configurations of the
form $\theta=\theta_{0}+2\pi nQ/L$ have topological charge $Q$ and action $0$,
as long as the system is large enough so that $\delta<2\pi Q/L$. When $\delta$
is sufficiently small, it is impossible to continuously222Here ”continuous”
means ”changing one spin at a time”. deform the action and change topological
charge without having an action equal to $s_{0}$, at least once. In [14] this
parameter was set to $s_{0}\rightarrow\infty$, creating an infinite barrier
between different topological sectors. While topological freezing usually
arises because of an _effective_ breaking of ergodicity, local algorithms like
Metropolis-Hastings and heat-bath break _explicitly_ ergodicity, when the
barrier is infinite. The present formulation of SDB is designed to accelerate
the mixing of local update algorithms and will not be able to deal with that
problem. As a consequence, we choose $s_{0}$ to be finite making barriers
between topological sectors potentially large, but not infinite.

##  3 Skewed Detailed Balance - SDB

###  3.1 Global Balance and Detailed Balance

Let us denote the probability of a given configuration $c$ by ${\cal
P}\left(c\right)$ and the transition probabilities, defining a given MCMC
algorithm, for the probability of transitioning from a configuration $c$ to
$c^{\prime}$, to be $T\left(c\rightarrow c^{\prime}\right)$. In order to
simulate correctly the model, we need to guarantee that the stationary
distribution of the Markov chain corresponds to our target probability
distribution,

${\cal P}\left(c\right)=\sum_{c^{\prime}}{\cal
P}\left(c^{\prime}\right)T\left(c^{\prime}\rightarrow c\right),$ (2)

which we shall refer as _global balance_. Proving that a given transition
matrix $T$ satisfies the above condition is a complicated task, in general. In
contrast, _detailed balance_ condition constitutes a more stringent
requirement

${\cal P}\left(c\right)T\left(c\rightarrow c^{\prime}\right)={\cal
P}\left(c^{\prime}\right)T\left(c^{\prime}\rightarrow c\right).$ (3)

Together with ergodicity, it guarantees that the probability distribution
simulated is precisely given by ${\cal P}$, which can be easily seen by
summing over $c^{\prime}$ on both sides of Eq. (3). A natural question is
whether we can use the freedom to break this condition to construct more
efficient algorithms.

###  3.2 Lifting and Skewed Detailed Balance

Skewed detailed balance condition provides an alternative path to satisfy Eq.
(2). The first step consists on applying the concept of lifting [6, 11, 15],
where we create a new replica of the phase space and introduce a _lifting
variable_ $\varepsilon=\\{+,-\\}$ that identifies the two copies. They occur
with the same probability $\varepsilon$: ${\cal
P}\left(c,\varepsilon\right)={\cal P}\left(c\right)/2$. In the enlarged space,
global balance reads ${\cal
P}\left(c\right)=\sum_{c^{\prime},\varepsilon^{\prime}}{\cal
P}\left(c^{\prime}\right)T\left(c^{\prime},\varepsilon^{\prime}\rightarrow
c,\varepsilon\right)$, We construct a transition matrix dependent on
$\varepsilon$ such that the _skewed balance_ condition is satisfied

${\cal P}\left(c\right)T\left(c,\varepsilon\rightarrow
c^{\prime},\varepsilon\right)={\cal
P}\left(c^{\prime}\right)T\left(c^{\prime},-\varepsilon\rightarrow
c,-\varepsilon\right).$ (4)

In order to make the algorithm ergodic in the enlarged space, the new variable
$\varepsilon$ must also be updated, i.e., there must be non-zero matrix
elements of the form $T\left(c,\varepsilon\rightarrow
c,-\varepsilon\right)\neq 0$. In contrast to the detailed balance condition,
Eq. (4) does not guarantee the correct probability distribution. By imposing
the latter and using Eq. (4) on the global balance condition, we arrive at the
condition

$T\left(c,\varepsilon\rightarrow
c,-\varepsilon\right)-T\left(c,-\varepsilon\rightarrow
c,\varepsilon\right)=\sum_{c^{\prime}\neq
c}\left(T\left(c,-\varepsilon\rightarrow
c^{\prime},-\varepsilon\right)-T\left(c,\varepsilon\rightarrow
c^{\prime},\varepsilon\right)\right),$ (5)

where we have implicitly assumed that $c$ and $\varepsilon$ are never updated
together. The equation above admits a large range of solutions. Here we will
always adopt the Turintsyn-Chertkov-Vucelja type solution [16], which is given
by

$T\left(c,\varepsilon\rightarrow
c,-\varepsilon\right)=\max\bigg{\\{}0,\sum_{c\neq
c^{\prime}}\left[T\left(c,-\varepsilon\rightarrow
c^{\prime},-\varepsilon\right)-T\left(c,\varepsilon\rightarrow
c^{\prime},\varepsilon\right)\right]\bigg{\\}}.$ (6)

Note that $\sum_{c\neq c^{\prime}}T\left(c,\varepsilon\rightarrow
c^{\prime}\varepsilon\right)$ is the probability of updating the configuration
$c$ to a different configuration $c^{\prime}$ at constant $\varepsilon$. With
this choice, as long as this probability is larger than its counterpart at
constant $-\varepsilon$, the system will remain in $\varepsilon$. This is a
nice property of the algorithm: it detects when the acceptance rate is larger
in the other replica, and only then there is a finite probability to move
there.

###  3.3 From Detailed Balance to Skewed Balance: The Skewness Function

Metropolis or heat-bath type of updates provide a straightforward way to
construct algorithms that satisfy detailed balance. A priory, it is not
obvious how to construct updates satisfying Eq. (4). Following [17], we can
adapt an algorithm satisfying detailed balance to construct an algorithm that
satisfies skewed balance on the extended model, after lifting has been
applied. In order to see how, assume that we have an algorithm satisfying
detailed balance, characterized by a transition matrix $T_{\mathrm{DB}}$ and a
_skewness function_ $\vartheta$ satisfying
$\vartheta_{\varepsilon}\left(c,c^{\prime}\right)=\vartheta_{-\varepsilon}\left(c^{\prime},c\right)$.
We can then verify that

$T\left(c,\varepsilon\rightarrow
c^{\prime},\varepsilon\right)=\vartheta_{\varepsilon}\left(c,c^{\prime}\right)T_{\mathrm{DB}}\left(c\rightarrow
c^{\prime}\right)$ (7)

satisfies the skewed balance condition Eq. (4). The task of constructing the
algorithm is now shifted to find a suitable skewness function. Naturally, for
$T$ to be a proper transition matrix, it should remain positive and smaller or
equal to one. A sufficient condition is to impose
$0\leq\vartheta_{\varepsilon}\left(c,c^{\prime}\right)\leq 1,\
\forall\varepsilon c,c^{\prime}$. The function $\vartheta$ works as a post
selection factor, which will favor certain updates over others, depending on
the value of $\varepsilon$. This can be constructed by taking into account the
effect of the update on an observable of choice [18, 19, 11]. Suppose we take
an observable $O$ that, for each configuration $c$,
$O\left(c\right)\in\mathbb{R}$. We can then define

$\vartheta_{\varepsilon}\left(c,c^{\prime}\right)=\frac{1+\varepsilon\lambda\mathrm{sign}\left(O\left(c^{\prime}\right)-O\left(c\right)\right)}{1+\lambda},$
(8)

where we have used the sign function, which returns the sign of its argument,
and $\lambda$ is a free parameter that can be used to explore the degree to
which we break detailed balance. If $\lambda=0$, the skewness function is
always 1, and detailed balance is satisfied. Conversely, if $\lambda=1$ all
updates that change the observable in a direction opposed to $\varepsilon$
will get rejected.

This approach has been used, for Ising and Potts models, with the
magnetization defining the skewness function [18, 19, 11]. The idea consists
of constructing an algorithm that can move from smaller to larger values of
magnetization by eliminating some of the randomness of the original one. If
$\lambda=1$, then we either only accept increasing the magnetization (when
$\varepsilon=1$) or decreasing it (when $\varepsilon=-1$) until the acceptance
rate has become small and the replica variable is updated. This has been shown
to provide an effective advantage.

##  4 Instanton-Skewed Detailed Balance - ISDB

As previously discussed, topological freezing occurs because different
topological sectors become separated by increasingly large barriers. Breaking
detailed balance may provide a pathway to overcome these barriers: if over a
succession of updates, we always move in a specific direction, we may hope to
climb a large barrier and update the topological sector. It is worth pointing
out how this scenario differs from the one describing magnetization in spin
models. The magnetization can be suitably updated locally, which is not the
case for the topological charge (except for possible lattice artifacts).
Furthermore, different values of magnetization are not separated by large
action barriers. The most obvious choice of observable for the skewness
function would be the topological charge $Q$ or the topological susceptibility
$Q^{2}/V$, where $V$ is the volume of the system. Because we will be using a
heat-bath algorithm, which is a local update algorithm, we need a different
choice for the skewness function.

###  4.1 Instantons and the Skewness Function

As previously discussed, for the $1-$d $O\left(2\right)$ model, one can define
configurations that minimize the action: $\theta_{n}^{Q}=nQ\frac{2\pi}{L}$. We
would like to construct a skewness function that promotes the transition
between different topological charges. We propose defining a family of
observables that measure how close the configuration is to a given
$\theta^{Q}$. We define

$O_{Q}\left(c\right)=\sum_{n}\cos\left(\theta_{n}^{Q}-\theta_{n}\right),$ (9)

where the values $\theta_{n}$ correspond to the angular values, in the
configuration $c$. We do not expect the algorithm to be crucially dependent on
the specific choice of the observable, as long as it constitutes a reasonable
measure of how far the two configurations differ from each other. This
specific choice is proportional to the (standard) action resulting from the
interaction of the system with an instanton in the background. We will
investigate the algorithm where the skewness function is defined as

$\vartheta_{\varepsilon}\left(c,c^{\prime}\right)=\frac{1+\varepsilon\lambda\mathrm{sign}\left(O_{1}\left(c^{\prime}\right)-O_{1}\left(c\right)\right)}{1+\lambda}.$
(10)

This will amount to promote the algorithm to approach the instanton
configuration given by $\theta_{n}=\theta_{n}^{Q=1}$ , if $\varepsilon=1$, and
$\theta_{n}=\pi+\theta_{n}^{Q=1}$ otherwise. We will call this algorithm
_instanton-skewed detailed balance_ algorithm. The goal of the remainder of
this paper is to take a standard heat-bath algorithm (equivalent to set
$\lambda=0$ in (10)), satisfying detailed balance, and compare it to its
instanton-skewed detailed balance adaptation with $\lambda=1$.

##  5 Updating with ISDB

The proposed algorithm biases the moves to copy one of two instantons on the
background, both with charge 1. This looks unnatural, and one may question if
the MCMC algorithm even converges to the Boltzmann distribution defined by the
action (1). As argued in the previous section, since SDB is satisfied (Eq.
(4)) and transitions between replicas are suitably chosen (Eq. (6)), this is
guaranteed to be the case. In this section, we discuss some of the
consequences of this algorithm and how it can constitute the building block of
more efficient algorithms. The results presented will focus on the constraint
angle $\delta=0.75$, system sizes $L\in\\{64,128\\}$ and barrier heights
$s_{0}\in\\{10,100,1000\\}$.

###  5.1 Updating the Topological Charge

In Section 2, we have argued that $s_{0}$ represents the height of the
barriers that separate topological sectors. The larger its value, the harder
will be to update the topological sector. This can be seen explicitly in Fig.
1, where we plot the value of the topological charge, after equilibrating, as
a function of Monte Carlo time, for different barrier heights $s_{0}$, and two
different volumes.

Figure 1: Sampling the topological charge, over Monte Carlo time, for the
algorithms that satisfy and break detailed balance. By increasing the action
barrier, $s_{0}$, updating the topological sector becomes less frequent.

The plot shows explicitly how increasing the height of the barrier makes it
less likely to update the topological charge, both for the standard heat-bath
algorithm ($\lambda=0$) and for its detailed balance-breaking version
($\lambda=1$). There is no stark distinction between the algorithms that are
exhibited by this plot. One may even expect that the ISDB algorithm will end
up performing worse, as it favors transitions between two specific
configurations. The main question that we are interested in is whether it can
accelerate relaxation when starting near those specific configurations. If so,
the principles of construction of this algorithm may be used to construct a
more efficient algorithm that does not rely on such specificities.

###  5.2 Equilibration of Charged Configurations

In this preliminary stage, we focus in a simple question: how does the system
relax, when topological charge is present on an initial configuration? To
address this question, we consider a family of initial configurations defined
by

$\theta_{n}^{\mathrm{init}}\left(\theta_{0}\right)=\theta_{0}+\theta_{n}^{Q},$
(11)

with $\theta_{0}\in\left[0,\pi\right]$. Given such an initial configuration,
we investigate how the expectation value of the topological charge $Q$ evolves
with Monte Carlo time $\tau_{\mathrm{MC}}$. Once the system has reached
equilibrium, the expectation value of the topological charge is zero. We
perform several independent runs to determine how the value
$\left<Q\right>\left(\tau_{\mathrm{MC}}\right)$ depends on the value of
$\theta_{0}$. The expectation is that, at least for $\theta_{0}=0$, the ISDB
algorithm should reach $\left<Q\right>\left(\tau_{\mathrm{MC}}\right)\simeq 0$
for smaller $\tau_{\mathrm{MC}}$ than the standard algorithm. This is indeed
what is observed in Fig. 2.

Figure 2: Comparison between the detailed balance breaking algorithm in shades
of blue ($\lambda=1$) and its detailed balance counterpart in red
($\lambda=0$) for $L=128$. The different blue lines represent different
initial configurations, where $\theta_{0}$ varies in steps of $0.1$. The
efficiency of the standard algorithm is independent of the initial condition
and is represented by a single line. On the left
$\theta_{0}\in\left[0,\pi/2\right]$ and on the right
$\theta_{0}\in\left[\pi/2,\pi\right]$.

In red we plot the evolution of the topological charge in the standard heat-
bath algorithm, which takes the initial topological charge given to the
system, and monotonically evolves it to zero. By construction of the
algorithm, it does not depend on $\theta_{0}$. This is slower than the ISDB
algorithm for an initial $\theta_{0}=0$, which is expected. The ISDB algorithm
performs updates with an instanton background corresponding precisely to
$\theta_{0}=0$ endowing it with an ability to quickly move on from that
configuration. The more we deviate from that configuration, by increasing
$\theta_{0}$, the more the algorithm gets slower. The performance becomes
worse than the standard algorithm, for large times, once $\theta_{0}\gtrapprox
0.2$, and keeps getting worse until $\theta_{0}=\pi/2$. By increasing more
$\theta_{0}$, the configuration continues to differ more and more from the
instanton background in the replica $\varepsilon=1$, but starts approaching
the background instanton of the $\varepsilon=-1$, and the performance improves
with it.

##  6 Conclusion

Topological freezing is a fundamental problem in lattice QCD, which extends to
other physical theories. This typically occurs on local updating algorithms,
which have difficulty climbing topological barriers. Solutions to this problem
for certain models are known. For example, cluster algorithms can solve this
problem for $O\left(N\right)$ models in $N-1$ dimensions [12, 13]. Proposals
to solve this problem, when good cluster algorithms are not known, include
change boundary conditions, master-field simulations, machine learning
techniques, proposal of independent non-local updates or parallel tempering
(some examples are [20, 21, 22, 23, 24]).

Here we propose an alternative direction by exploring algorithms that break
detailed balance. While this is being formulated as an alternative to the
aforementioned methods, they are perfectly compatible. We have explored a
particular instance of breaking detailed balance, skewed detailed balance,
properly modified to improve sampling of the topological charge in the 1-d
$O\left(2\right)$ model, which we named Instanton-Skewed Detailed Balance
algorithm. This algorithm operates by keeping instantons "in the background".
By analyzing how the system relaxes from a configuration with topological
charge, we observe that if the initial configuration is close enough to one of
the possible instantons on the background, the system relaxes faster. The
observed behavior suggests that the proposed algorithm, in its current form,
is globally less efficient than the standard heat bath. However, for certain
configurations that are similar to the background instantons, the converse is
true. This gives hope that, by increasing the number of possible background
configurations, one might be able to outperform the original algorithm
altogether.

We emphasize that this is just a very specific way of breaking detailed
balance. The results presented here suggest that skewed detailed balance can
accelerate the sampling of different topological sections, but more research
on it and other extensions are necessary to clarify a potential role in
properly tackling the topological freezing problem.

## References

  * [1] B. Alles, G. Boyd, M. D’Elia, A. Di Giacomo and E. Vicari, _Hybrid monte carlo and topological modes of full qcd_ , _Physics Letters B_ 389 (1996) 107.
  * [2] L. Del Debbio, G.M. Manca and E. Vicari, _Critical slowing down of topological modes_ , _Physics Letters B_ 594 (2004) 315.
  * [3] S. Schaefer, R. Sommer, F. Virotta, A. Collaboration et al., _Critical slowing down and error analysis in lattice qcd simulations_ , _Nuclear Physics B_ 845 (2011) 93.
  * [4] G. McGlynn and R.D. Mawhinney, _Diffusion of topological charge in lattice qcd simulations_ , _Physical Review D_ 90 (2014) 074502.
  * [5] C. Bonati and M. D’Elia, _Topological critical slowing down: variations on a toy model_ , _Physical Review E_ 98 (2018) 013308.
  * [6] P. Diaconis, S. Holmes and R.M. Neal, _Analysis of a nonreversible markov chain sampler_ , _Annals of Applied Probability_ (2000) 726.
  * [7] E.P. Bernard, W. Krauth and D.B. Wilson, _Event-chain monte carlo algorithms for hard-sphere systems_ , _Physical Review E_ 80 (2009) 056704.
  * [8] H. Suwa and S. Todo, _Markov chain monte carlo method without detailed balance_ , _Physical review letters_ 105 (2010) 120603.
  * [9] R.D. Schram and G.T. Barkema, _Monte carlo methods beyond detailed balance_ , _Physica A: Statistical Mechanics and its Applications_ 418 (2015) 88.
  * [10] M. Kaiser, R.L. Jack and J. Zimmer, _Acceleration of convergence to equilibrium in markov chains by breaking detailed balance_ , _Journal of statistical physics_ 168 (2017) 259.
  * [11] F. Faizi, G. Deligiannidis and E. Rosta, _Efficient irreversible monte carlo samplers_ , _Journal of Chemical Theory and Computation_ 16 (2020) 2124.
  * [12] U. Wolff, _Collective monte carlo updating for spin systems_ , _Physical Review Letters_ 62 (1989) 361.
  * [13] U. Wolff, _Asymptotic freedom and mass generation in the o (3) nonlinear $\sigma$-model_, _Nuclear Physics B_ 334 (1990) 581.
  * [14] W. Bietenholz, U. Gerber, M. Pepe and U.-J. Wiese, _Topological lattice actions_ , _Journal of High Energy Physics_ 2010 (2010) 1.
  * [15] F. Chen, L. Lovász and I. Pak, _Lifting markov chains to speed up mixing_ , in _Proceedings of the thirty-first annual ACM symposium on Theory of computing_ , pp. 275–281, 1999.
  * [16] K.S. Turitsyn, M. Chertkov and M. Vucelja, _Irreversible monte carlo algorithms for efficient sampling_ , _Physica D: Nonlinear Phenomena_ 240 (2011) 410.
  * [17] Y. Sakai and K. Hukushima, _Eigenvalue analysis of an irreversible random walk with skew detailed balance conditions_ , _Physical Review E_ 93 (2016) 043318.
  * [18] Y. Sakai and K. Hukushima, _Dynamics of one-dimensional ising model without detailed balance condition_ , _Journal of the Physical Society of Japan_ 82 (2013) 064003.
  * [19] K. Hukushima and Y. Sakai, _An irreversible markov-chain monte carlo method with skew detailed balance conditions_ , in _Journal of Physics: Conference Series_ , vol. 473, p. 012012, IOP Publishing, 2013.
  * [20] M. Lüscher and S. Schaefer, _Lattice qcd without topology barriers_ , _Journal of High Energy Physics_ 2011 (2011) 1.
  * [21] M. Lüscher, _Stochastic locality and master-field simulations of very large lattices_ , in _EPJ Web of Conferences_ , vol. 175, p. 01002, EDP Sciences, 2018.
  * [22] G. Kanwar, M.S. Albergo, D. Boyda, K. Cranmer, D.C. Hackett, S. Racaniere et al., _Equivariant flow-based sampling for lattice gauge theory_ , _Physical Review Letters_ 125 (2020) 121601.
  * [23] D. Albandea, P. Hernández, A. Ramos and F. Romero-López, _Topological sampling through windings_ , _The European Physical Journal C_ 81 (2021) 1.
  * [24] C. Bonanno, M. D’Elia, B. Lucini and D. Vadacchino, _Towards glueball masses of large-n su(n) yang-mills theories without topological freezing via parallel tempering on boundary conditions_ , _arXiv preprint arXiv:2210.07622_ (2022) .


# Cost models for geo-distributed massively parallel streaming analytics

Anna-Valentini Michailidou and Anastasios Gounaris
Aristotle University of Thessaloniki, Greece
<EMAIL_ADDRESS>Konstantinos Tsichlas
University of Patras, Greece
<EMAIL_ADDRESS>
###### Abstract

This report is part of the DataflowOpt project on optimization of modern
dataflows and aims to introduce a data quality-aware cost model that covers
the following aspects in combination: (1) heterogeneity in compute nodes, (2)
geo-distribution, (3) massive parallelism, (4) complex DAGs and (5) streaming
applications. Such a cost model can be then leveraged to devise cost-based
optimization solutions that deal with task placement and operator
configuration.

##  1 Introduction

A fundamental research activity with a view to devising new methodologies
regarding advanced dataflow (i.e., data-intensive workflow) optimization
solutions is to develop and adopt appropriate cost models. Such cost models
need to be capable of accurately reflecting the response time in parallel
environments, taking into account issues such as partial overlapping, CPU core
sharing and resource contention that are all common in modern dataflow
executions. Furthermore, the cost models need to be extended to cover
additional dataflow objectives, such as data quality. Finally, these cost
models need to cover arbitrarily complex DAG (directed acyclic graph)
execution plans. This report, apart from examining cost models that are
directly amenable to optimizations, further discusses tractability issues
regarding the optimization problems that typically accompany cost models.

More specifically, in Section 2, we present alternative approaches to develop
cost models with a view to covering the landscape that is characterized by the
combination of (1) heterogeneity in compute nodes, (2) geo-distribution, (3)
massive parallelism, (4) complex DAGs and (5) streaming applications. We then
highlight the need to develop our own cost model, which covers the afore-
mentioned aspects, i.e., heterogeneity, massive parallelism and geo-
distribution, while being suitable for data-intensive streaming applications.
This cost model, is presented in detail in Section 3 and the relevant code is
in https://github.com/annavalentina/Equality.222This code repository contains
broader work. Section 4 concludes the report.

##  2 Cost Models for Dataflow Optimization

In this section, we discuss the main representatives regarding the cost models
for data intensive analytics, along with the associated optimization problems,
whenever such problems accompany a cost model presentation, and tractability
analysis. The aspects that are of immediate interest include: massive
parallelism, heterogeneity in terms of resources and (geo-)distribution, while
we assume that the complete job is modelled by a (potentially complex) DAG.
Please note that we use the terms DAG, topology, workflow (or dataflow)
interchangeably in the text.

###  2.1 A Cost Model for a Multi-core Shared Memory Setting Tailored for
Streaming Data Analytics

Zhang et al. [37] propose a cost model that is used to maximize the processing
throughput of data stream processing systems on shared-memory multi-core
architectures. The throughput $R$ is equal to the sum of the output rates of
all the sink nodes of the topology, that is $R=\sum_{sink}{r_{o}}$. The output
rate $r_{o}$ for a time interval $t$ is estimated using the number of tuples
that are processed, $num$, and the time needed to process them $t_{p}$. Hence
$r_{o}=num/t_{p}$. $num$ is equal to the number of tuples generated from the
previous nodes in the topology, $num=\sum_{producers}{r_{i}*t}$ where $r_{i}$
is the input rate. The $t_{p}$ variable is equal to
$t_{p}=\sum_{producers}{r_{i}*t*T(p)}$ where $T(p)$ is the average time spent
on handling each tuple and consists of the fetching ($T^{f}$) and execution
($T^{e}$) times. The fetching time $T^{f}$ is 0 in case the producer keeps its
data locally; otherwise it is equal to $[N/S]*L_{(i,j)}$ where $N$ is the
average size of a tuple, $S$ the cache line size and $L_{i,j}$ the worst case
memory access latency between two CPU sockets $i$ and $j$ in a generic NUMA
computing architecture. The variables $L_{i,j}$ and $S$ are given as
statistical input metadata, whereas $N$ and $T^{f}$ are found using profiling
techniques. Finally, the input rate of the source operators is also given by
the problem setting (denoted as $I$). The algorithms proposed in [37] aim to
maximize the throughput $R$ by deciding a placement of nodes to CPU sockets
and a replication level for each operator. The proposed algorithms focus on
the degree of parallelism of the operators taking into account the
heterogeneity of the resources, without considering geo-distribution. Apart
from this difference from our main focus, the heterogeneity in [37] is limited
and refers to the locality or not of data in a NUMA architecture, which
impacts on the fetching time as described above.

In summary, this cost model emphasizes on the different fetching times based
on the data memory address and the placement of the operator.

####  2.1.1 Optimization Problem Details and Tractability

Zhang et al. [37] aim to maximize the application processing throughput, given
an input stream rate by placing each operator on specific CPU sockets $S_{i}$,
where $i=1\dots n$, and deciding an optimal replication level in a multicore
shared-memory architecture. The formulation of the problem they try to solve
is the following:

$maximize\;\sum_{sink}\overline{r_{o}}(i.e.,~{}expected~{}r_{o})$
$s.t.,\forall i,j\\\ \ \epsilon\\\ \ 1,...,n,$ $\sum_{operators\\\ \ at\\\ \
S_{i}}\overline{r_{o}}*T\leq C,$ (1) $\sum_{operators\\\ \ at\\\ \
S_{i}}\overline{r_{o}}*M\leq B,$ (2) $\sum_{operators\\\ \ at\\\ \
S_{j}}\sum_{producers\\\ \ at\\\ \ S_{i}}\overline{r_{o}(s)}*N\leq Q_{i,j}.$
(3)

Equation (1) ensures that the aggregated CPU demand on each socket will not
exceed its CPU capacity $C$, where $T$ is the average time spent on handling
each tuple. Similarly, equation (2) ensures that the aggregated amount of
bandwidth requested on a socket (calculated using the average memory bandwidth
consumption per tuple, $M$) will not exceed its maximum attainable local DRAM
bandwidth $B$. Finally, equation (3) enforces that the data transfer from one
socket to another must not exceed the maximum remote channel bandwidth. It is
also constrained that each operator is placed exactly once.

Each operator of the total $|o|$ operators can be replicated at most $k$
times, thus the total possible configurations are $k^{|o|}$. For each of these
configurations, there exist $m^{n}$ different placements, where $m$ is the
number of sockets and $n$ the number of replicas. Due to the exponential
search space and the corresponding ILP that needs to be solved, a brute-force
approach is prohibitive. The authors propose a branch and bound algorithm that
iteratively places operators under a given replication level and then tries to
increase the degree of replication for the bottleneck operator. In total,
there exist $(n*m)^{n}$ candidate solutions ($n*m$ for $n$ iterations) and
three heuristics are used in order to speed up the traversal of the solution
space.

###  2.2 Cost Models for Data Flows Tailored to Parallel Executions

Kougka et al. [20] present three cost models for data flow response time in a
massively parallel environment taking into account execution overlaps. Their
approach can be extended to parallel distributed settings but heterogeneity of
the resources is not considered.

The first two models refer to pipelined segments where the tasks in the data
flow DAG form a chain and, also, no movement of data over the network occurs,
since the execution takes place on a single machine (thus making communication
cost time negligible). The idea behind the first cost model is that each task
can be executed on a separate core resulting to overlapping execution times.
In such a case, the slowest execution time defines the total response time.
Also the authors in [20] account for overheads coming from factors such as
memory locks, contention and context switching between threads, where all
these factors are captured by a single variable $\alpha$. The model proposed
is the following: $ResponseTime(RT)=\alpha*max\\{c_{1},...,c_{n}\\}$, where
$c_{i}$ is the execution time of task $i$ in the DAG.

In the second model, the case where a core needs to take on more than one task
is examined. The new model is the following:
$ResponseTime(RT)=\alpha*max\\{max\\{c_{1},...,c_{n}\\},\sum\\{c_{1},...,c_{n}\\}/m\\}$
where $m$ is the total number of cores. The outcome of this model is the same
as the first one in case there exists a single task with execution time higher
than the sum of the execution times of all the tasks.

Finally, a generalized model is presented that takes into account multiple
pipeline segments and multiple machines. The model presented is the following:
$ResponseTime(RT)=\sum z_{i}*w^{c}*c_{i}+\sum z_{ij}*w^{cc}*cc_{i\rightarrow
j}$. The first sum refers to execution time whereas the second one to
communication cost time. The $z$ variables are binary ones that are equal to 1
when the task and the communication costs contribute to the final $RT$,
capturing the overlapping of tasks. The $w$ variables generalize the $\alpha$
variable used in the previous models and capture the possible overheads of the
execution. Finally, $cc_{i\rightarrow j}$ represents the communication cost of
the edges between the tasks.

In summary, the cost models in [20] emphasize on capturing the effect of
concurrent task executions on the response time in a parallel homogeneous
setting.

In a similar setting, Ali et al. [4] propose a cost model for optimizing User
Defined Functions (UDFs) in data-intensive workflows, such as ETL (Extract-
Transform-Load) workflows. They try to optimize the execution of the UDFs in
ETL workflows based on the execution time and monetary costs of task
invocations. The optimization is divided into three stages. In the first one,
the cost model determines if a UDF can be parallelized and in the second
stage, it determines the degree at which it should be parallelized. Finally,
in the third stage, the cost model generates an efficient machine
configuration. The cost model is used to simulate the UDF execution and derive
estimates regarding execution time as well as other metrics like the dataset
size, the execution speed, CPU, I/O, memory characteristics etc.

In summary, the cost model in [4] focuses on the impact of increased degree of
parallelism regarding certain expensive operators in the DAG.

####  2.2.1 Optimization Problem Details and Tractability

The cost model in [20] mentioned above is leveraged to find an optimal
(partial) ordering of the operators (nodes) in a DAG, in the case that
multiple orderings can produce the same result. The response time optimization
in such a setting is not simpler than the general task ordering optimization
problem aiming at the minimization of the sum of the costs, as reported in
[21]. This problem, as proven in [8], is intractable and more importantly “it
is unlikely that any polynomial time algorithm can approximate the optimal
plan to within a factor of $O(n^{\theta})$”, where $\theta$ is some positive
constant and $n$ is the number of the tasks in the data flow.

Ali et al. [4], as already explained, optimize the execution of UDFs in ETL
workflows using also the sum cost metric. The workflows are divided into
stages and each stage can be processed using multiple code/implementation
variants. The objective is to find the optimal variant for each stage. In
total, there exist $n^{m}$ possible variants ($m$ stages, $n$ code variants
for each stage). This problem is NP-hard since the Multiple Choice Knapsack
Problem (MCKP) is a special case of it. In the MCKP, given $m$ stages
$N_{1},...,N_{m}$ of the workflow to pack in a knapsack of capacity $B$ (the
budget constraint), where each code variant $j\;\epsilon\;N_{i},i=1,2,...,m$,
has execution time $T_{i,j}$ and cost $C_{i,j}$, the problem is to choose
exactly one code variant from each stage such that the sum of execution times
is minimized while the sum of costs does not exceed the budget $B$. The
proposed formulation is the following (the decision variable is $x_{ij}$):

$minimize\;(Z)\qquad\sum_{i=1}^{m}\sum_{j\epsilon N_{i}}T_{i,j}*x_{i,j}$

Subject to:

$\sum_{i=1}^{m}\sum_{j\epsilon N_{i}}C_{i,j}*x_{i,j}\leq B,$ (1)
$\sum_{j\epsilon N_{i}}x_{i,j}=1,i=1,...,m,$ (2)
$x_{i,j}\;\epsilon\;\\{0,1\\},j\;\epsilon\;N_{i},i=1,...,m.$ (3)

###  2.3 Cost Models Tailored for Streaming Data Analytics in an Edge/
Fog/Cloud Computing Setting

Hiessl et al. [15] aim to minimize response time, migration and enactment cost
as well as to maximize availability in data stream processing systems by
placing (and replacing) operators of stream topologies. However, the operator
execution is not parallelized, i.e., each DAG node is placed on exactly one
compute node. Their proposal targets fog environments, where heterogeneity and
geo-distribution of the resources are inherent characteristics and they extend
the system model and formulation presented in [9]. The adopted cost model
represents both the resources and operators as graphs. In the former case, the
links represent the network links where, in the latter case, the graphs refer
to the flow of data streams between the operators. Each compute node has its
CPU frequency, number of CPU cores, memory and storage capacity, whereas each
operator has its own CPU frequency, number of CPU cores, memory and storage
requirements. Also, for each operator, there exists a subset of the resources
that are available to allocate its execution. The operators can be enacted and
migrate from one compute node to another. The response time of a topology is
equal to the maximum delay of all paths in the stream data flow, which depends
on the processing time of tuples and on the delay of sending the tuples. The
modeling includes the $C_{op}(x)$ variable, which refers to the total
enactment cost per second, which is the cost of running the topologies on the
available resources. Another variable is $C_{mig}(x)$, which is the cost that
comes from all the migrations of operators from one compute node to another
and depends on the operator image size as well as the data rate of pulling
this operator to a compute node.

In general, this model focuses on capturing the environment heterogeneity and
then, in the subsequent optimization problem, the emphasis is placed on
expressing several constraints.

Renart et al. [29] also place operators across edge and cloud devices in order
to optimize metrics like end-to-end latency, WAN traffic and messaging cost.
The setting they work on is geo-distributed with heterogeneous resources but
the operators are not parallelized, as in the cost model above. Each resource
is characterized by its CPU and memory capabilities ($cpu,mem$). Each network
link $l_{k\leftrightarrow l}$ has a certain bandwidth capability
$bdw_{k\leftrightarrow l}$ with other resources. Each operator has its own CPU
and memory requirements ($cpu,mem$), selectivity, process rate $\mu_{<i,k>}$,
input and output size $\varsigma^{in}_{i},\varsigma^{out}_{i}$ and output
production rate $\lambda^{out}_{i}$, which is based on its input rate
$\lambda^{in}_{i}$. The resource and operator modeling is based on an M/M/1
queue model, which allows for calculating the computation and communication
times. $stime_{<i,k>}=\frac{1}{\mu_{<i,k>}-\lambda^{in}_{i}}$ gives the
computation time of an operator $i$ placed on a resource $k$ and
$ctime_{<i,k><j,l>}=\frac{1}{(\frac{bdw_{k\leftrightarrow
l}}{\varsigma^{out}_{i}})-\lambda^{in}_{j}}+l_{k\leftrightarrow l}$ denotes
the communication time for an operator $i$ placed on resource $k$ to send data
to operator $j$ placed on resource $l$. The end-to-end latency $L_{p_{i}}$ of
a path $p_{i}$ is the sum of the computation times of all the operators in
this path plus the communication times between them. The way the latency of
each path is calculated in [29] is similar to that of [15] but the latter
focuses on optimizing the maximum latency of all paths, whereas this work
focuses on the sum of the latencies. The WAN traffic of a path $W_{p_{i}}$
indicates the size of messages that cross it. The messaging cost of a path
$C_{p_{i}}$ shows the number of messages between edge and cloud. Finally the
total aggregate cost model they use is the following:
$AggregateCost_{p_{i}}=w_{l}*\frac{L_{p_{i}}}{Par_{lat}}+w_{w}*\frac{W_{p_{i}}}{Par_{wan}}+w_{c}*\frac{C_{p_{i}}}{Par_{cost}}$,
where $w_{l}+w_{w}+w_{c}=1$ are weights and $Par_{lat},Par_{wan},Par_{cost}$
indicate the current latency, WAN traffic and messaging cost of the running
application.

This cost model is characterized by the usage of queue theory and the fact
that it aims to capture the environment heterogeneity, which is a main goal
for cost models in geo-distributed settings. As previously, the final cost is
a weighted sum of multiple individual costs, which are computed by simple
summation expressions, given that these cost models are typically leveraged by
multi-objective placement optimization solutions.

####  2.3.1 Optimization Problem Details

Hiessl et al. [15] place stream operators in fog computing environments. Both
the computing resources and the application are presented using topology
graphs. The aim is to map the operators and edges from the application graph
to the resources and network links of the resource graph, respectively, in
order to minimize the following objective function:
$F^{\prime}_{cost}=F^{\prime}(x,y,r)+w_{c_{op}}\frac{C_{opmax}-C_{op}(x)}{C_{opmax}-C_{opmin}}+w_{c_{mig}}\frac{C_{migmax}-C_{mig}(x)}{C_{migmax}-C_{migmin}}$,
where
$F^{\prime}(x,y,r)=w_{r}\frac{R_{max}-r}{R_{max}-R_{min}}+w_{a}\frac{logA(x,y)-logA_{min}}{logA_{max}-logA_{min}}$.
The $w$ variables are weights used in the simple additive weighting method
that applies to multi-objective problems. The $F^{\prime}(x,y,r)$ variable
focuses on response time and availability ($R$ and $A$ respectively).

The formulation of the problem is the following:

$minimize\;F^{\prime}_{cost}$

Subject to:

$C_{op}(x)\leq B_{op}$ (1) $C_{mig}(x)\leq B_{mig}$ (2) $\sum_{u\epsilon
V_{res}}T_{(actual,i)}\frac{S_{i}}{S_{u}}x_{i,u}\leq T_{(max,i)}\qquad\forall
i\;\epsilon\;operators$ (3) $\sum_{i\epsilon
V_{dsp}}P_{(CPU,i)}P_{(Cores,i)}x_{i,u}\leq
P_{(CPU,u)}P_{(Cores,u)}\qquad\forall u\;\epsilon\;resources$ (4)
$\sum_{i\epsilon V_{dsp}}P_{(Mem,i)}x_{i,u}\leq P_{(Mem,u)}\qquad\forall
u\;\epsilon\;resources$ (5) $\sum_{i\epsilon V_{dsp}}P_{(HD,i)}x_{i,u}\leq
P_{(HD,u)}\qquad\forall u\;\epsilon\;resources$ (6)
$r\geq\sum_{k=1}^{n_{p}}\sum_{u\epsilon
V^{i_{k}}}\frac{ET_{i}}{S_{u}}x_{i_{k},u}+\sum_{k=1}^{n_{p}-1}\sum_{(u,v)\epsilon\;available\\_resources}d_{(u,v)}y_{(i_{k},i_{k+1}),(u,v)}\qquad\forall
p\;\epsilon\;paths$ (7)

Equations (1) and (2) limit the enactment and migration cost using budget
constraints $B_{op}$ and $B_{mig}$, respectively. Equation (3) ensures that
the processing duration of an operator assigned to a resource should not
exceed the maximum processing time limit; in the equation, $S_{i}$ denotes the
speed of processing of operator $i$ in previous optimization circles and
$S_{u}$ the speed of processing of resource $u$. Equations (4), (5) and (6)
ensure that the CPU, Memory and storage requirements of the operators are
fulfilled while not exceeding the capacities of the resources. Finally,
equation (7) ensures that the total response time is greater than or equal to
the response times of all the topology paths; $d$ denotes the network delay
and $ET$ the execution time per tuple, To avoid quadratic or higher order
equations, the authors skip the multiplication of the decision variables
$x_{i,u}$ and $y_{(i,j),(u,v)}$ which refer to the placement of operators and
links, respectively.

Renart et al [29] place operators of IoT dataflows on either edge or cloud
resources in order to minimize end-to-end latency, data transfer and messaging
cost between the edge and the cloud. The aggregate cost formula for a single
path has been presented above. The formulation of the problem they try to
solve that is given below refers to the aggregation costs across all paths:

$minimize\;\sum_{p_{i}\epsilon P}AggregateCost_{p_{i}}$

Subject to:

$\lambda^{in}_{o}<\mu_{<o,r>}\qquad\forall o\;\epsilon\;operators,\forall
r\;\epsilon\;resources|mo_{<o,r>}=1$ (1)
$\lambda^{in}_{o}<(\frac{bdw_{k\leftrightarrow
n}}{\varsigma^{out}_{o-1}})\qquad\forall o\;\epsilon\;operators,\forall
k\leftrightarrow n\;\epsilon\;links|mo_{<o,k>}=1$ (2) $\sum_{o\epsilon
O}mo_{<o,r>}*\lambda^{in}_{o}\leq cpu_{r}\qquad\forall r\;\epsilon\;resources$
(3) $\sum_{o\epsilon O}mo_{<o,r>}*mem_{<o,r>}\leq mem_{r}\qquad\forall
r\;\epsilon\;resources$ (4) $\sum_{\begin{subarray}{c}s_{i}\rightarrow\epsilon
S\\\ k\leftrightarrow l\epsilon L\end{subarray}}ms_{<i\rightarrow
j,k\leftrightarrow l>}*\varsigma^{out}_{i}\leq bdw_{k\leftrightarrow
l}\qquad\forall k\leftrightarrow l\;\epsilon\;links$ (5) $\sum_{r\epsilon
R}mo_{<o,r>}=1\qquad\forall o\;\epsilon\;operators$ (6)
$\sum_{k\leftrightarrow l\epsilon L}ms_{<i\rightarrow j,k\leftrightarrow
l>}=1\qquad\forall s_{i\rightarrow j}\;\epsilon\;streams$ (7)

The variable $mo_{<o,r>}$ is set to 1 if operator $o$ is placed on resource
$r$ and, similarly, $ms_{<i\rightarrow j,k\leftrightarrow l>}=1$ if the stream
between operators $i,j$ is assigned on the link between resources $k,l$.
Equations (1) and (2) ensure that the input rate does not exceed the process
rate of the resource and that the link is not saturated, respectively.
Equations (3) and (4) ensure that the CPU and memory requirements of each
operator are fulfilled and the resource capacities are not exceeded. Equation
(5) ensures that the bandwidth limits of the links are not exceeded. Finally,
equations (6) and (7) ensure that an operator is placed on only one resource
and that each data stream is placed on only one link, respectively.

####  2.3.2 Tractability

In both [15] and [29], the problem is more complex than the problem of
choosing among multiple execution engines, when the placement decision of a
node affects the other placements. This problem is analyzed in [22], and is
proven that it is NP-hard, and moreover, _cannot be approximated by a
polynomial-time algorithm with approximation error bound less than 8/7_. Also,
evidence about the NP-hardness is provided in [9]. As shown from the problem
formulations, such task placement problems are specific forms of constrained
optimizations (mixed) integer linear programming problems.

###  2.4 A Cost Model Tailored to Stage-by-stage Analytics in a Cloud
Computing Setting

Gounaris et al. [13] propose a bi-objective cost model in order to optimize
queries in multi-cloud environments. The queries are in a form of DAGs and are
divided into strides where each stride represents a step in the execution and
consists either of one operator or a group of operators running in parallel.
However each operator cannot be parallelized. The goal is to assign strides to
cloud-based VMs. Each cloud provider offers VMs to users and charges based on
one of the following policies : on-demand, where the payment is computed based
on the exact usage time, reserved, where the user pays upfront for a specific
time and then pays on-demand, spot instance, which is a bidding process or
committed, where the customer rents the VM for a predefined time. The model
they present estimates the execution time and monetary cost and models the
charging policies and fees. The basic formula for modeling the execution time
is the following: $TotalTime=\sum_{s=1}^{n}max(S_{s,1\rightarrow
i}^{VM_{k\rightarrow k^{\prime}}},...,S_{s,m^{s}\rightarrow
j}^{VM_{k\rightarrow k^{\prime}}})$, where $n$ is the number of strides,
$m^{s}$ the number of operators in the $s^{th}$ stride and $S_{s,i\rightarrow
j}^{VM_{k\rightarrow k^{\prime}}}=O_{s,i}^{VM_{k}}+T_{s,i\rightarrow
j}^{VM_{k\rightarrow k^{\prime}}}$. $O_{s,i}^{VM_{k}}$ is the time it takes to
execute the $i^{th}$ operator of the $s^{th}$ stride on $VM_{k}$ and depends
on the hardware characteristics of the heterogeneous VMs. $T_{s,i\rightarrow
j}^{VM_{k\rightarrow k^{\prime}}}$ is the time to transfer data produced from
operator $i$ which runs on $VM_{k}$ to operator $j$ on $VM_{k^{\prime}}$. In
case there exist network bottlenecks, the total time is not the max of the sum
of execution and transmission time of each operator but rather the sum of them
as follows: $TotalTime=\sum_{s=1}^{n}\sum_{i=1}^{m^{s}}S_{s,i\rightarrow
j}^{VM_{k\rightarrow k^{\prime}}}$. In case pipelining is supported, an
operator can start transmitting data before its execution is complete. Then
the next formulation holds $S_{s,i\rightarrow j}^{VM_{k\rightarrow
k^{\prime}}}=max(O_{s,i}^{VM_{k}},T_{s,i\rightarrow j}^{VM_{k\rightarrow
k^{\prime}}}$). The monetary cost of a query is the following:
$\sum_{s=1}^{n}\sum_{i=1}^{m^{s}}Price(S_{s,i\rightarrow j}^{VM_{k\rightarrow
k^{\prime}}},F_{pr}(P_{t}^{\alpha},VM_{k}),F_{pr}(P_{t^{\prime}}^{\alpha},VM_{k^{\prime}}))$,
where $F_{pr}$ represents the price offer from cloud provider $pr$ for
$Vm_{k}$ according to the pricing policy $P_{t}^{\alpha}$. Then Price finds
the fee for using $VM_{k}$ for $S$ time and transferring data to
$VM_{k^{\prime}}$.

####  2.4.1 Optimization Problem Details and Tractability

This cost model accompanies bi-objective problems initially encountered in
federated databases, where an execution plan is processed stage-by-stage
across multiple sites, which are shown to be NP-hard [28, 34].

###  2.5 A Cost Model for (Incremental) MapReduce Jobs

Li et al. [23] propose a system built on top of incremental Hadoop that
allocates resources and uses runtime scheduling techniques in order to meet
user latency requirements while maximizing throughput. The queries are
represented using DAGs where each node is either a data distributor or a Map
Reduce pair. They model the mean and the variance of latency both per-tuple
and per-window. Regarding the per-tuple latency, G/G/1 queues are used to
represent streams. The total latency and variance are a sum of latencies and
variances due to 12 causes including batching, queuing, heartbeats awaiting,
CPU, network and Disk I/O. The mean CPU latency per tuple is formulated as
follows:
$E(\ddot{L}_{cpu})=\frac{u}{2*min(1-p,1/n)*C}$, where $u$ is the total
resource required by a batch and depends on the average number of tuples and
the average resource required per tuple. Also, $p$ is the fraction of the
resource required by the other threads on the same node since resources
receive multiple map-reduce tasks, $n$ is the number of cores of a CPU and $C$
is the number of processing cycles of a CPU per unit time. Network and Disk
I/O are modeled in a similar way. The queuing delay can be calculated based on
the theory of the G/G/1 queues. The per-window latency is modeled in cases
where the results of a window are computed only when a reducer has received
all the tuples and is equal to $E(\tilde{L})=E(U)+E(F)$ where $U$ is the
maximum of all the tuples’ latencies in this window and F represents the
execution of a particular partitioned window.

In summary, this cost model focuses on splitting the latency factors of a
tuple processed by a MapReduce round into 12 parts, and takes into account
resource sharing.

Work |  Setting |  Objective |  Resources’ metadata required |  Operators’ metadata required
---|---|---|---|---
Zhang et al. [37] |  Multicore shared-memory architectures: decide placement of operators to sockets and degree of parallelism |  Processing Throughput |  CPU cycles per socket, DRAM bandwidth, bandwidth between sockets, memory access latency, cache line size |  Memory bandwidth consumption, execution time, average size of tuples, input stream rate
Kougka et al. [20] |  Distributed multicore machines: decide operator ordering |  Response Time |  Communication cost, number of cores |  Execution cost, selectivity
Hiessl et al. [15] |  Edge/Fog environment: decide placement of operators on compute nodes |  Response Time, availability, migration and enactment cost |  Availability, delay of network links, speedup, CPU frequency, number of CPU cores, memory and storage capacity per compute node |  Speedup from previous run, image size, requirements for CPU frequency, number of CPU cores, memory and storage per operator
Renart et al. [29] |  Edge/cloud environment: decide placement of operators on compute nodes |  End-to-end latency, WAN traffic, messaging cost |  CPU and memory capacities, bandwidth costs |  CPU and memory requirements, selectivity, input/output event rate, input/output event size
Gounaris et al. [13] |  Multi-cloud environments: decide placement of operators on compute nodes |  Execution time, monetary cost |  Charging policy, CPU, I/O access speed |  Size of data
Ali et al. [4] |  Multicore machines: decide degree of parallelism |  Execution time, monetary costs |  Configuration, row and byte process rate |  Size of data, memory, I/O and CPU usage
Li et al. [23] |  Multicore machines running MR: decide degree of parallelism (amount of resources) and scheduling |  Latency, throughput |  CPU, disk I/O and network characteristics |  Estimated data rate, resource required and input-to-output ration
Table 1: Summary of techniques examined.

####  2.5.1 Optimization Problem Details and Tractability

The work in Li et al. [23] decides the amount of resources per Map/Reduce
solving a non-linear constrained optimization problem (using MinConNLP). In
the generic case, these problems are intractable. The main focus is on
scheduling the tuples through re-ordering their execution on-the-fly.

###  2.6 Summary of Techniques

We summarize the techniques examined thus far in Table 1. We also provide an
overview of their main limitations, which call for the development of a new
cost model, tailored to massively parallel complex streaming analytics jobs
over heterogeneous geo-distributed nodes; this cost model is presented in
Section 3. The cost models in [37, 20, 4] do not cover heterogeneity in
compute node characteristics and/or locations. The models in [15, 29, 13] do
not account for massive parallelism. Therefore, they are inadequate and
unsuitable to be adopted for our case, which focuses on complex DAGs that are
executed over geo-distributed heterogeneous resources while benefiting from
partitioned parallelism.

###  2.7 Additional Issues

We have already identified that the existing cost models fall short in
supporting the use cases targeted hereby. For completeness, in this section,
we examine more cost modelling approaches and we also consider scheduling
problems.

####  2.7.1 Other Modeling Approaches

Karnagel et al. [19] estimate the performance of tasks in order to achieve an
optimal adaptive placement of OLAP queries on heterogeneous hardware. To
estimate the query runtime, they use an automated black-box online-learning
approach. They use monitoring for execution times as well as input data
cardinalities. To estimate the transfer costs, they rely on statistics
collected in the calibration phase rather than performing monitoring. In
general, learning task costs using (machine)-learning approaches is quite
common and is orthogonal to how individual costs are combined to estimate the
cost of the complete DAG.

Learning methods can also apply to the whole DAG. More specifically, Witt et
al. [36] survey predictive performance modeling (PPM) techniques that estimate
metrics like execution time, required memory and predict future wait times of
tasks and jobs in data-parallel and distributing systems. They compare the
works by the performance factors they consider; that is workload (resource
usage patterns), heterogeneity (of resources), scale (amount of resources and
problem size) and contention (caused by resource sharing). They study both
task performance models and job performance ones. The representative works we
present from this survey here consider the performance factors of
heterogeneity and scalability. In the first category (namely, task performance
models), a single task runs on a single node. Marin et al. [24] predict the
performance of tasks, the cost of the computation and the execution duration
using memory access patterns that are monitored using histograms, the cache
size and input size. Matsunaga et al. [25] predict execution duration, output
file size and resident set size of a task using regression trees and
characteristics like CPU speed, cache size, available memory and disk
capabilities, number of threads as well as location and size of input. Iverson
et al. [17] predict the execution duration of a task using the input size and
hardware characteristics and a local learning approach. The second category
refers to approaches that deal with multiple nodes that comprise a job.
Delimitrou et al. [11] predict the performance by observing other applications
with similar resource usage patterns and make use of Singular Value
Decomposition to classify incoming jobs. Shanjay et al. [30] use regression to
model the computation duration and the non-overlapped communication duration,
whereas the total execution duration is equal to their sum. Cunha et al. [10]
aim to predict performance in order to choose between cloud and a local
cluster. To achieve this, they propose a linear model and an empirical one
that favors local execution when the prediction is unclear. Last, Alipourfard
et al. [5] propose an off-line search method that uses a Gaussian model to
select cloud configurations like number of virtual machines, number of cores,
CPU speed, RAM per core etc in order to optimize the performance of a job.

Finally, the work of Agrawal et al. [2] selects the task ordering and the
assignment of tasks to compute nodes, while it targets both the bottleneck
cost metric (called period) and the critical path cost (which is equivalent to
the latency); this work considers precedence constraints between the tasks in
the DAG but the exact ordering is flexible. Each task has a cost and a
selectivity (selectivities are independent). Also, there is a set of servers
acting as the compute nodes. The focus is on both homogeneous and
heterogeneous platforms (heterogeneity refers to network speeds) but not on
geo-distributed settings. The communication speed is considered in some
flavors. If there are no communication costs, for the cases of homogeneous
machines, Agrawal et al. present algorithms that perform optimal ordering in
polynomial time. They also prove that if the machines are heterogeneous, the
problems are NP-hard even when no precedence constraints exist. On the other
hand, if there are communication costs, they consider the models of
overlapping (communication fully overlaps with computation), and sequential
(no overlap) processing. The problem is always NP-complete, apart from the
case of period minimization and overlap. In general, in this work the cost
models are simple and the value lies in the algorithmic analysis and
solutions, which however refer to a more restricted setting than the one
assumed in this report.

####  2.7.2 Cost Models in Scheduling Works

We organize the discussion in this section around five recent surveys in the
area of scheduling.

First, Hong et al. [16] review algorithms, infrastructures and architectures
for resource management on edge/fog computing. Infrastructure refers to
hardware, software and middleware. The architectures are divided into data
flow, control and tenancy ones and they present data flow architectures that
provide techniques which take into consideration bandwidth, latency, energy
efficiency, quality, security and heterogeneity in order to reduce the
communication cost. Regarding the algorithms, they divide them into four
categories; discovery, benchmarking, load-balancing and placement. Discovery
algorithms aim to identify edge resources to use for deployment of
applications. Benchmarking focuses on capturing the actual performance (CPU,
memory, storage) using stress tests and tools. Load-balancing aims to
distribute tasks using several approaches, such as optimization techniques,
graph-based balancing, breadth-first search or cooperative load balancing.
Finally, placement algorithms assign tasks onto resources and can be divided
into dynamic, condition-aware and iterative ones.

The placement techniques in [16] are summarized as follows. Wang et al. [35]
propose off-line and on-line techniques for finding an optimal placement of
task (service) instances on Mobile Micro Clouds in order to minimize the
average cost while predicting future cost parameters with known accuracy. The
off-line algorithm finds the optimal placement of services within a predefined
time-window whereas the on-line approximation one makes placement decisions
upon arrival of the services. Taneja et al. [33] map IoT application jobs on
fog resources by ordering both based on their CPU, RAM and bandwidth
requirements and capacities, respectively. In case, it is not feasible to map
jobs to fog nodes, these are delegated to cloud resources. Souza et al. [32]
propose three algorithms, namely First-Fit, Best-Fit and Best-Fit with Queue,
for placing services, which are divided into smaller atomic-services that form
DAGs, on fog or cloud nodes. First-fit assigns services to the first available
fog node, Best-fit first orders the services based on their requirements and
then assigns them. Both these algorithms assign services to the cloud, in case
no available fog node exists, whereas Best-fit with Queue puts these services
in queues while waiting for fog nodes to become available by predicting if
that will be more beneficial. Fernandoet al. [12] propose Honeybee, a model
for load balancing in a group of heterogeneous mobile nodes by using the work-
stealing method. The independent jobs are pulled by the mobile sources instead
of being pushed to them and, additionally, a node can “steal" jobs from
another node. Also, this work supports the mobility of the devices, in the
sense that they can join or leave the group at any time. Skarlat et al. [31]
place IoT services on fog nodes taking into account their QoS requirements
(e.g. execution deadlines). They organize the fog devices, called fog cells,
into fog colonies, which act as data centers. Each colony has a control node
used for management and communication with the other colonies and the cloud.
The placement problem is formulated as an Integer Linear Programming one. Ni
et al. [27] propose a dynamic resource allocation strategy for fog computing
based on Priced Time Petri Nets (PTPN) in order to improve the resource
utilization while satisfying the users’ QoS requirements. The authors also
propose algorithms to predict the time cost and price cost of the tasks as
well as compute the credibility evaluation of resources. The users can select
the resources from a set of provided resources based on these characteristics.
All these proposals do not account for multiple compute nodes co-operating for
processing the same task in complex DAGs.

In a second survey, Beaumont et al. [7] present and compare algorithms that
focus on scheduling tasks on CPU and GPU resources. The scheduling algorithms
decide which task to execute first, when to start the execution and where to
allocate it. They provide performance guarantees in the form of lower makespan
bounds. They categorize the algorithms as on-line and off-line and the tasks
as independent or with precedence constraints. In the on-line scenario they
prove lower bound conditions, while, in the off-line algorithms, a linear
programming problem is solved to keep the workload on CPU and GPU processors
under a defined threshold and also to satisfy any precedence relations between
tasks. Nevertheless, these solutions cannot be transferred to massively
parallel (streaming) settings, while also heterogeneity of inter-node
communication speeds are not considered. On the other hand, online scheduling
techniques can inspire adaptive resource placement policies in the context of
our project.

Figure 1: Workflow orchestration challenges taxonomy proposed in [6] .

The third survey considered is authored by Barika et al. [6], who review the
challenges in orchestrating big data analysis workflows. First, they present
the different programming models that exist, namely MapReduce, Stream, NoSQL,
SQL, Message Ingestion and Hybrid. They also provide a taxonomy for the
challenges, which can be divided into three categories: cloud-based, data-
based and workflow-based ones as depicted in Figure 1. Next, they focus on
some of the current approaches and techniques for the workflow-based
challenges. In a nutshell, according to their discussion, the workflow
specification language can be generic or custom with the latter limiting the
ease to use multiple execution environments. The initialization of workflow
execution aims to partition the workflow into fragments. This can be done by
taking into consideration constraints like security, privacy, storage and
compute capacities, data transfer constraints, or a combination of them or
taking into consideration only the task and data dependencies within the
workflow. Parallelization can be coarse-grained focusing on workflow level
(when dealing with multiple workflows) or fine-grained focusing on
activity/task level. Scheduling can be push-based, which aims to schedule
tasks to available workers or pull-based, where the workers request for the
tasks they will process. The former can be further categorized into static,
dynamic or hybrid. Workflow fault-tolerance can be achieved either through
reactive or proactive techniques. Finally the security challenge can be
tackled by using multi-cloud architectures or replication-level techniques. No
cost models are explicitly discussed, but the significance to cover issues
such as resource, storage and data transmission characteristics is mentioned.

Adhikari et al. [1] survey scheduling strategies for workflows in cloud
environments. They offer equations and definitions for multiple optimization
objectives. The execution time of a task is equal to
$ET_{ki}=\frac{SZ_{k}}{CP_{i}}$ where $SZ_{k}$ is the size of the task and
$CP_{i}$ the capacity of the cloud VM instance which is the product of the
number of cores and the size of each core (in MIPS or MFLOPS). The completion
time of a task is either equal to its execution time or equal to the sum of
its execution and waiting times. The waiting time for a task is 0 or equal to
the maximum completion time of its predecessor tasks. The execution cost of a
task refers to the cost of running on a VM instance and is equal to
$\frac{ET_{ki}}{\tau}*CO_{i}$ where $\tau$ is the time period the task used
the VM and $CO_{i}$ represents the cost of the VM. The energy consumption of a
task refers to the amount of energy required to run on a VM instance and is
equal to $E_{ik}=DP_{ik}+SP_{ik}$, that is the sum of the dynamic and static
power consumption. The makespan (total execution time) of a workflow is the
maximum completion time of its tasks, the total execution cost is the sum of
the execution costs of all the tasks and similarly the total energy
consumption is the sum of the energy consumption of all the tasks. Some of the
most common scheduling objectives are total execution time (makespan), total
execution cost, schedule length ratio (SLR-the ratio of the completion time of
a task to the total execution time of the workflow), availability,
reliability, energy consumption, resource utilization, communication overhead
and security. The authors also discuss some of the emerging trends in workflow
scheduling which include containers (e.g. Docker, Kubernetes), server-less
computing (DynamoDB, Amazon SNS, Google cloud pub/sub) and fog computing. The
latter environment also raises the need to consider the latency objective
between two compute nodes $i$ and $j$ ,which is equal to
$LT_{i,j}=PD_{i,j}+SD_{i,j}$. $PD_{i,j}=\frac{DS_{i,j}}{BW_{i,j}}$ is the
propagation delay and depends on the distance between the nodes $DS_{i,j}$ and
the bandwidth of the network $BW_{i,j}$. $SD_{i,j}=\frac{SZ_{i,j}}{TR_{i,j}}$
is the serialization delay and depends on the size of the workflow in bits and
the transmission rate of the network. However, the content of this survey is
not suitable for our setting, because it does not account for placement for
complex streaming applications, where all DAG nodes are active simultaneously
and partitioned parallelism is employed.

Finally, Salaht et al. [3] survey techniques for the service placement problem
in Fog Computing environments. They characterize the works as being
centralized or decentralized, off-line or on-line, static or dynamic and as to
whether they support mobility of the devices or not. Furthermore, they focus
on the optimization problems (multi- vs. mono-objective), the metrics
considered (latency, resource utilization, cost, energy consumption,
congestion ratio, etc.), the formulation of the problem (integer linear
programming, integer non-linear programming, mixed integer linear programming,
mixed integer quadratic programming) and finally the resolution strategies
(exact solution, heuristics, meta-heuristics, approximations). In terms of
modeling, the constraints to be covered include resource, network and
applications (e.g., locality of a task placement) ones. The main limitation of
the techniques examined that account for complex DAGs to be applicable to our
setting is that they do not consider massive parallelism.

##  3 The Proposed Cost Model Proposal

The cost model we propose is tailored to streaming data scenarios focusing on
geo-distributed, heterogeneous resources. The heterogeneity apart from network
speeds also refers to hardware characteristics, like CPU and RAM. Also we
consider the massive parallelism of the operators. Since it refers to a
streaming case, in terms of performance, the focus is on latency. The cost
model extends the one presented in [14] which however, apart from smaller
differences, does not cover massive parallelism on top of arbitrary resources.

Symbol |  Meaning
---|---
$G_{op}$ |  Graph representing the overall analytics job
$V_{op}$ |  Vertices of $G_{op}$ (operators or tasks)
$E_{op}$ |  Edges of $G_{op}$ (data shuffling)
$ED$ |  Edge devices
$s_{i}$ |  Selectivity of $i\in V_{op}$
$available_{i,u}$ |  availability of $u\in ED$ for $i\in V_{op}$
$comCost_{u,v}$ |  Communication cost between $u,v\in ED$
$ED_{i}\subset ED$ |  Subset of edge devices $i\in V_{op}$ can be assigned to
$x_{i,u}$ |  Fraction of $i\in V_{op}$ assigned on $u\in ED$
$\alpha$ |  Network congestion factor
$enabledLinks_{i,j}$ |  Number of devices that participate in $(i,j)\in E_{op}$
Table 2: Notations used in the proposed cost model.

Table 2 presents the main notation used. Specifically, we represent a
streaming analytics job as a Directed Acyclic Graph (DAG)
$G_{op}=(V_{op},E_{op})$ where each node in $V_{op}$ represents an operator
and the edges in $E_{op}$ represent the data flow between operators. An
operator represents a set of analysis job steps at a fine level of granularity
that can run on the same device (typically in a pipelined manner). Each
operator can be physically partitioned into multiple instances, where each
instance is responsible for a disjoint data partition. Edges represent data
re-distribution among the partitioned operator instances. For example,
according to the above description, in Apache Spark, operators correspond to
Spark job stages while edges are placed when data shuffling takes place. The
data sources are the operators in the $G_{op}$ without incoming edges, and we
assume that they produce data in batches periodically.

Each operator takes as input multiple tuples (i.e., records) in batches that
need to be analyzed and also has a certain selectivity $s_{i}$ according to
its functionality. For example, a transformation operator with selectivity
$s_{i}=1$ outputs a tuple for each tuple it takes as input, while a filtering
operator with $s_{i}=0.5$ outputs a tuple for every two input tuples on
average. If operator $i$ is a sink node in the graph, then, there is no
selectivity, whereas, if it is a source node, its selectivity is equal to 1.
Each pair of edge devices has a different communication cost, $comCost_{u,v}$
($u,v\in ED$), which relates to the connected devices network speed; the
latter is also dependent on the physical distance and link capacity between
the devices.

Each operator is assigned to multiple edge devices that run in parallel. That
means that each device $u$ is assigned a fraction $x_{i,u}$ of tuples of an
operator $i$ to analyze. However, due to privacy, security and other reasons,
there exists a subset of edge devices where an operator can be assigned to
(denoted by the flag $available_{i,u}$ for operator $i$ regarding device $u$).
We define as $ED{i}\subset ED$ the subset of edge devices operator $i$ can be
assigned to.

As already stated, the QoS metric we try to minimize is the average latency.
This latency is equal to the latency of the critical path (i.e., the slowest
path) with regards to a single input data batch and consists of the average
communication latency between the operators in the critical path. In order for
the critical path’s latency to be representative of the whole graph’s latency,
we assume there is no waiting time between the tasks of the operators. This is
achieved through setting the operator requirements at such levels, so that
they sustain the input data production rate. A natural consequence of this is
that each operator can start its execution as soon as it receives its input
batch of tuples from its parent nodes in the graph. As our focus is on geo-
distributed realms, we make the realistic assumption that the execution
latency of each operator is negligible and the communication cost dominates.
The total latency of a tuple is the time to flow from its source node
downstream to a sink node.

The communication latency between two nodes is expressed by
$max\\{x_{i,u}*s_{i}*\sum_{v\in ED_{j}}\\\ comCost_{u,v}*x_{j,v})\\},~{}u\in
ED_{i}$ across all instances of operator $i$. This is equal to the slowest
data transfer that comes from a single device and refers to the cost of a
batch of input data. However, this modeling is over-simplistic because it does
not take into account the overhead of an operator instance maintaining
multiple remote connections. In order to take into account such an overhead,
we introduce a parameter, notated as $\alpha$ multiplied by the number of
enabled links. A link between device $u$ and device $v$ for two operators
$i,j$ is enabled when $x_{i,u}\neq 0$, $x_{j,v}\neq 0$ and $u\neq v$. Thus,
$enabledLinks_{i,j}$ denotes the number of devices that exchange data between
two operators over the network. The total latency of the topology is equal to
the latency of the critical path, that is slowest path of the graph that leads
from a source node to a sink one (not including). We denote as $path$ any path
from a source to an operator just upstream a sink operator ($j$ is the
operator just after $i$).

$Latency=max_{path\in G_{op}}\\{\sum_{i\in V_{op}\in
path}max\\{x_{i,u}*s_{i}*\\\ \sum_{v\in
ED_{j}}(comCost_{u,v}*x_{j,v})+\alpha*enabledLinks_{i,j}\\}\\}$

###  3.1 Accounting for Additional Objectives

The above cost model targets a single objective. It is trivial, through simple
sum functions to consider the movement of data over the network, as in [26],
or the total time resources are occupied for processing a specific time
period.

An important novelty, however, in modern dataflows if their focus on
additional objectives, such as quality optimization. The quality of the data
is an important aspect in IoT scenarios. Low quality can lead to less useful
results. The quality of data can be categorized into multiple dimensions. Some
examples of these are _completeness, timeliness_ and _accuracy_. Some of the
most common factors that lead to decreased data quality include the
heterogeneity of data sources, missing and dirty data due to network
malfunctions or security constraints [18]. Consider the scenario where a user
wants to analyze the output of a sensor every 5 minutes. If, for some reason,
the sensor malfunctions or restarts, the output will be inaccurate or
misleading and the user should ignore it. Thus there is a need to “rate" the
quality of the data and/or to detect outliers and missing values that could
lead to misleading results. However, this data quality rating can be very time
consuming especially when dealing with streaming big data that need to be
processed under tight latency constraints.

$DQ_{fraction}$ denotes the percentage of the input data on which the quality
check is applied. Then, we can combine the quality and latency metrics in a
single formula as follows:

$F=\frac{Latency}{1+\beta*(DQ_{fraction})},~{}\beta\geq 0$ (8)

The latency represents the average time a tuple needs to come out of a sink
node, beginning from a source one. The latency is increased with more data
quality checks but the relationship is not proportional. This is due to the
fact that the more the quality checks, the less an edge device can be assigned
tasks of upstream operators, thus inducing higher communication cost. The
configuration parameter $\beta$ is a weight denoting the importance of data
quality (DQ). The higher its value, the more beneficial (in terms of
minimizing $F$) becomes to increase $DQ_{fraction}$. Setting $\beta$ to 0
essentially removes DQ from the optimization criteria.

_Example._ Lets assume we have a simple linear DAG with 3 nodes (operators)
and 3 devices. The selectivity of each node is the following:
$s_{0}=1,s_{1}=1.5$ ($s_{2}$ does not have any impact). The communication
costs between devices are presented in Table 3. For simplicity we set $\alpha$
equal to 0.

Device | 0 | 1 | 2
---|---|---|---
0 | 0 | 1.5 | 2
1 | 1.5 | 0 | 1
2 | 2 | 1 | 0
Table 3: Communication cost between edge devices in GBps. operator/device | 0 | 1 | 2
---|---|---|---
0 | 0.8 | 0.2 | 0.0
1 | 0.7 | 0.0 | 0.3
2 | 0.3 | 0.4 | 0.3
Table 4: Fraction of operator assigned to each edge device.

We assign the operators as shown in Table 4. The transfer time for each link
in the graph depends on the selectivity of that operator, the communication
cost of the devices and the assigned fractions. Each device sends data to
every other device based on the fraction of the current and the succeeding
operator. For example, for operator 0, device 0 will take $0.8\times 1\times
0\times 0.7$ sec/unit to send data to device 0 (itself), $0.8\times 1\times
1.5\times 0$ sec/unit to device 1 and $0.8\times 1\times 2\times 0.3$ sec/unit
to device 2. We sum these cost elements to derive the total cost of
communication between operator 0 and operator 1 for the first device with id
0; this sum is 0.48. For device 1, we have $0.2\times 1\times(1.5\times
0.7+0\times 0+1\times 0.3)$ = 0.27 sec/unit. For device 2, the communication
cost is 0. Therefore, the latency due to the link $0\rightarrow 1$ is
$max\\{0.48,0.27,0\\}=0.48$. Similarly, we can derive that the cost between
$1\rightarrow 2$ is $max\\{1.26,0,0.45\\}=1.26$ sec/unit. Therefore, the
overall latency is 1.74 sec/unit.

If we assume that the $DQ_{fraction}$ is 0.5 then with $\beta=1$, the
objective function $F$ becomes 1.16. Further assume that we examine another
scenario, where $DQ_{fraction}$ is 1 at the expense of moving the complete
fraction $x_{2,0}$ to device 2, i.e., the last operator runs 40% on device 1
and 60% on device 2. Then, the latency cost of $1\rightarrow 2$ becomes
$max\\{1.89,0,0.18\\}=1.89$ and the complete latency is 2.37. The new value of
$F$ is 1.185, i.e., despite the important increase in latency, the new plan is
not better in terms of $F$. However, if we give even more weight to the data
quality through increasing the $\beta$ value, e.g., setting $\beta=2$, then
the initial and modified allocation yield $F$ values of 0.87 and 0.79,
respectively, i.e., the trade-off in the modification has become beneficial.

##  4 Summary

In this report, we first describe existing cost modeling approaches
highlighting their inadequacy to cover all aspects relevant to a modern
setting including heterogeneity in compute node characteristics and network
connections, massive parallelism, and complex DAG-shaped streaming application
execution plans. In addition, we briefly discuss tractability issues related
to the optimization problems associated with cost model proposals. To cover
the aforementioned gap, in Section 3, we propose a new cost model, which also
accounts for data quality.

_Acknowledgment._ The research work was supported by the Hellenic Foundation
for Research and Innovation (H.F.R.I.) under the “First Call for H.F.R.I.
Research Projects to support Faculty members and Researchers and the
procurement of high-cost research equipment grant” (Project Number:1052,
Project Name: DataflowOpt).

## References

  * [1] Mainak Adhikari, Tarachand Amgoth, and Satish Narayana Srirama.  A survey on scheduling strategies for workflows in cloud environment and emerging trends.  ACM Comput. Surv., 52(4):68:1–68:36, 2019.
  * [2] Kunal Agrawal, Anne Benoit, Fanny Dufossé, and Yves Robert.  Mapping filtering streaming applications.  Algorithmica, 62(1-2):258–308, 2012.
  * [3] Farah Aït-Salaht, Frédéric Desprez, and Adrien Lebre.  An overview of service placement problem in fog and edge computing.  ACM Comput. Surv., 53(3):65:1–65:35, 2020.
  * [4] Syed Muhammad Fawad Ali and Robert Wrembel.  Towards a cost model to optimize user-defined functions in an ETL workflow based on user-defined performance metrics.  In Tatjana Welzer, Johann Eder, Vili Podgorelec, and Aida Kamisalic Latific, editors, Advances in Databases and Information Systems - 23rd European Conference, ADBIS 2019, Bled, Slovenia, September 8-11, 2019, Proceedings, volume 11695 of Lecture Notes in Computer Science, pages 441–456. Springer, 2019.
  * [5] Omid Alipourfard, Hongqiang Harry Liu, Jianshu Chen, Shivaram Venkataraman, Minlan Yu, and Ming Zhang.  Cherrypick: Adaptively unearthing the best cloud configurations for big data analytics.  In Aditya Akella and Jon Howell, editors, 14th USENIX Symposium on Networked Systems Design and Implementation, NSDI 2017, Boston, MA, USA, March 27-29, 2017, pages 469–482. USENIX Association, 2017\.
  * [6] Mutaz Barika, Saurabh Kumar Garg, Albert Y. Zomaya, Lizhe Wang, Aad van Moorsel, and Rajiv Ranjan.  Orchestrating big data analysis workflows in the cloud: Research challenges, survey, and future directions.  ACM Comput. Surv., 52(5):95:1–95:41, 2019.
  * [7] Olivier Beaumont, Louis-Claude Canon, Lionel Eyraud-Dubois, Giorgio Lucarelli, Loris Marchal, Clément Mommessin, Bertrand Simon, and Denis Trystram.  Scheduling on two types of resources: A survey.  ACM Comput. Surv., 53(3):56:1–56:36, 2020.
  * [8] Jen Burge, Kamesh Munagala, and Utkarsh Srivastava.  Ordering pipelined query operators with precedence constraints.  Technical Report 2005-40, Stanford InfoLab, 2005.
  * [9] Valeria Cardellini, Vincenzo Grassi, Francesco Lo Presti, and Matteo Nardelli.  Optimal operator placement for distributed stream processing applications.  In Avigdor Gal, Matthias Weidlich, Vana Kalogeraki, and Nalini Venkasubramanian, editors, Proceedings of the 10th ACM International Conference on Distributed and Event-based Systems, DEBS ’16, Irvine, CA, USA, June 20 - 24, 2016, pages 69–80. ACM, 2016.
  * [10] Renato Luiz de Freitas Cunha, Eduardo R. Rodrigues, Leonardo P. Tizzei, and Marco A. S. Netto.  Job placement advisor based on turnaround predictions for HPC hybrid clouds.  Future Gener. Comput. Syst., 67:35–46, 2017.
  * [11] Christina Delimitrou and Christos Kozyrakis.  Quasar: resource-efficient and qos-aware cluster management.  In Rajeev Balasubramonian, Al Davis, and Sarita V. Adve, editors, Architectural Support for Programming Languages and Operating Systems, ASPLOS ’14, Salt Lake City, UT, USA, March 1-5, 2014, pages 127–144. ACM, 2014.
  * [12] Niroshinie Fernando, Seng W. Loke, and Wenny Rahayu.  Computing with nearby mobile devices: A work sharing algorithm for mobile edge-clouds.  IEEE Trans. Cloud Comput., 7(2):329–343, 2019.
  * [13] Anastasios Gounaris, Zisis Karampaglis, Athanasios Naskos, and Yannis Manolopoulos.  A bi-objective cost model for optimizing database queries in a multi-cloud environment.  J. Innov. Digit. Ecosyst., 1(1-2):12–25, 2014.
  * [14] Lin Gu, Deze Zeng, Song Guo, Yong Xiang, and Jiankun Hu.  A general communication cost optimization framework for big data stream processing in geo-distributed data centers.  IEEE Trans. Computers, 65(1):19–29, 2016.
  * [15] Thomas Hiessl, Vasileios Karagiannis, Christoph Hochreiner, Stefan Schulte, and Matteo Nardelli.  Optimal placement of stream processing operators in the fog.  In 3rd IEEE International Conference on Fog and Edge Computing, ICFEC 2019, Larnaca, Cyprus, May 14-17, 2019, pages 1–10. IEEE, 2019.
  * [16] Cheol-Ho Hong and Blesson Varghese.  Resource management in fog/edge computing: A survey on architectures, infrastructure, and algorithms.  ACM Comput. Surv., 52(5):97:1–97:37, 2019.
  * [17] Michael A. Iverson, Füsun Özgüner, and Lee C. Potter.  Statistical prediction of task execution times through analytic benchmarking for scheduling in a heterogeneous environment.  IEEE Trans. Computers, 48(12):1374–1379, 1999.
  * [18] Aimad Karkouch, Hajar Mousannif, Hassan Al Moatassime, and Thomas Noel.  Data quality in internet of things: A state-of-the-art survey.  Journal of Network and Computer Applications, 73:57 – 81, 2016\.
  * [19] Tomas Karnagel, Dirk Habich, and Wolfgang Lehner.  Adaptive work placement for query processing on heterogeneous computing resources.  Proc. VLDB Endow., 10(7):733–744, 2017.
  * [20] Georgia Kougka and Anastasios Gounaris.  Optimization of data flow execution in a parallel environment.  Distributed Parallel Databases, 37(3):385–410, 2019.
  * [21] Georgia Kougka, Anastasios Gounaris, and Alkis Simitsis.  The many faces of data-centric workflow optimization: a survey.  Int. J. Data Sci. Anal., 6(2):81–107, 2018.
  * [22] Georgia Kougka, Anastasios Gounaris, and Kostas Tsichlas.  Practical algorithms for execution engine selection in data flows.  Future Gener. Comput. Syst., 45:133–148, 2015.
  * [23] Boduo Li, Yanlei Diao, and Prashant J. Shenoy.  Supporting scalable analytics with latency constraints.  Proc. VLDB Endow., 8(11):1166–1177, 2015.
  * [24] Gabriel Marin and John M. Mellor-Crummey.  Cross-architecture performance predictions for scientific applications using parameterized models.  In Edward G. Coffman Jr., Zhen Liu, and Arif Merchant, editors, Proceedings of the International Conference on Measurements and Modeling of Computer Systems, SIGMETRICS 2004, June 10-14, 2004, New York, NY, USA, pages 2–13. ACM, 2004.
  * [25] Andréa M. Matsunaga and José A. B. Fortes.  On the use of machine learning to predict the time and resources consumed by applications.  In 10th IEEE/ACM International Conference on Cluster, Cloud and Grid Computing, CCGrid 2010, 17-20 May 2010, Melbourne, Victoria, Australia, pages 495–504. IEEE Computer Society, 2010.
  * [26] Anna-Valentini Michailidou and Anastasios Gounaris.  Bi-objective traffic optimization in geo-distributed data flows.  Big Data Res., 16:36–48, 2019.
  * [27] Lina Ni, Jinquan Zhang, Changjun Jiang, Chungang Yan, and Kan Yu.  Resource allocation strategy in fog computing based on priced timed petri nets.  IEEE Internet Things J., 4(5):1216–1228, 2017.
  * [28] Christos H. Papadimitriou and Mihalis Yannakakis.  Multiobjective query optimization.  In Proceedings of the Twentieth ACM SIGACT-SIGMOD-SIGART Symposium on Principles of Database Systems, May 21-23, 2001, Santa Barbara, California, USA. ACM, 2001.
  * [29] Eduard Gibert Renart, Alexandre Da Silva Veith, Daniel Balouek-Thomert, Marcos Dias de Assunção, Laurent Lefèvre, and Manish Parashar.  Distributed operator placement for iot data analytics across edge and cloud resources.  In 19th IEEE/ACM International Symposium on Cluster, Cloud and Grid Computing, CCGRID 2019, Larnaca, Cyprus, May 14-17, 2019, pages 459–468. IEEE, 2019.
  * [30] H. A. Sanjay and Sathish S. Vadhiyar.  Performance modeling of parallel applications for grid scheduling.  J. Parallel Distributed Comput., 68(8):1135–1145, 2008.
  * [31] Olena Skarlat, Matteo Nardelli, Stefan Schulte, and Schahram Dustdar.  Towards qos-aware fog service placement.  In 1st IEEE International Conference on Fog and Edge Computing, ICFEC 2017, Madrid, Spain, May 14-15, 2017, pages 89–96. IEEE Computer Society, 2017.
  * [32] Vitor Barbosa C. Souza, Xavier Masip-Bruin, Eva Marín-Tordera, Sergio Sánchez-López, Jordi Garcia, Guang-Jie Ren, Admela Jukan, and Ana Juan Ferrer.  Towards a proper service placement in combined fog-to-cloud (F2C) architectures.  Future Gener. Comput. Syst., 87:1–15, 2018.
  * [33] Mohit Taneja and Alan Davy.  Resource aware placement of iot application modules in fog-cloud computing paradigm.  In 2017 IFIP/IEEE Symposium on Integrated Network and Service Management (IM), Lisbon, Portugal, May 8-12, 2017, pages 1222–1228. IEEE, 2017\.
  * [34] Efthymia Tsamoura, Anastasios Gounaris, and Kostas Tsichlas.  Multi-objective optimization of data flows in a multi-cloud environment.  In Proceedings of the Second Workshop on Data Analytics in the cloud, pages 6–10, 2013.
  * [35] Shiqiang Wang, Rahul Urgaonkar, Ting He, Kevin Chan, Murtaza Zafer, and Kin K. Leung.  Dynamic service placement for mobile micro-clouds with predicted future costs.  IEEE Trans. Parallel Distributed Syst., 28(4):1002–1016, 2017\.
  * [36] Carl Witt, Marc Bux, Wladislaw Gusew, and Ulf Leser.  Predictive performance modeling for distributed computing using black-box monitoring and machine learning.  CoRR, abs/1805.11877, 2018.
  * [37] Shuhao Zhang, Jiong He, Amelie Chi Zhou, and Bingsheng He.  Briskstream: Scaling data stream processing on shared-memory multicore architectures.  In Peter A. Boncz, Stefan Manegold, Anastasia Ailamaki, Amol Deshpande, and Tim Kraska, editors, Proceedings of the 2019 International Conference on Management of Data, SIGMOD Conference 2019, Amsterdam, The Netherlands, June 30 - July 5, 2019, pages 705–722. ACM, 2019\.