# On the positive mass theorem in general relativity and Lorentz covariance
of the Dirac wave equation in quantum mechanics

Changbiao Wang<EMAIL_ADDRESS>ShangGang Group, 10 Dover Road, North
Haven, CT 06473, USA

###### Abstract

The positive mass theorem in general relativity states that in an
asymptotically flat spacetime, if the momentum–energy tensor is divergence-
free and satisfies a dominant energy condition, then a total momentum–energy
four-vector can be formed, of which the energy component is nonnegative. In
this paper, we take the wave four-tensor of a plane light wave in free space
as a counterexample to show that there is no guarantee that a total four-
vector can be formed. Thus the theoretical framework for the positive mass
theorem is flawed. In addition, it is also shown as well that the Lorentz
covariance of Dirac wave equation is not compatible with Einstein mass–energy
equivalence.

###### pacs:

03.50.De, 04.20.-q, 04.20.Cv

##  I Introduction

In the theory of relativity, the momentum and energy of a closed physical
system can be described by a four-tensor $T^{\mu\nu}$, usually called
momentum–energy tensor. It is well established in traditional textbooks (r1, ,
p. 443) (r2, , p. 46) (r3, , p. 166) (r4, , p. 756) that if the tensor is
divergence-free, namely $T^{\mu\nu}_{~{}~{}\,,\nu}=0$, then the total momentum
and energy of the system can be obtained by carrying out integration of the
time-column elements of the tensor to form a constant four-vector
$P^{\mu}=\int T^{\mu 4}\textrm{d}^{3}x$, which is usually called _conservation
law_. However recently, it has been shown by enumerating specific
counterexamples that this conclusion is not correct in general r5 ; r6 ; r7 ;
r8 .

In an asymptotically flat spacetime, if the momentum–energy tensor
$T^{\mu\nu}$ satisfies the conservation law $T^{\mu\nu}_{~{}~{}\,,\nu}=0$, and
also satisfies an additional _dominant energy condition_ , namely the energy
component $T^{44}\geq|T^{\mu\nu}|$ r9 , then a total momentum–energy four-
vector $P^{\mu}=\int T^{\mu 4}\textrm{d}^{3}x$ can be formed, of which the
total energy component $P^{4}=\int T^{44}\textrm{d}^{3}x$ is nonnegative,
which is called _positive mass theorem_ in the general relativity r10 ; r11 ;
r12 ; r13 ; r14 . (Note that Rindler’s practice of the tensor indices
$\mu,\nu=1,2,3,4$ (r15, , p. 138) is used throughout the paper except for in
Sec. II.)

The best way to disprove a mathematical conjecture is to give its
counterexample. In this paper, we take the wave four-tensor
$T^{\mu\nu}=K^{\mu}K^{\nu}$ as a counterexample to show that there is no
guarantee that a total four-vector can be formed although it satisfies the
traditional conservation law $T^{\mu\nu}_{~{}~{}\,,\nu}=0$ and the dominant
energy condition $T^{44}\geq|T^{\mu\nu}|$, where $K^{\mu}$ is the wave four-
vector of a plane light wave in free space, set up by Einstein r16 . From this
we conclude that the theoretical framework for the positive mass theorem is
fundamentally flawed.

Finally, it is also shown as well that the Lorentz covariance of Dirac wave
equation for an electron is not compatible with Einstein mass–energy
equivalence.

##  II Positive mass theorem

In this section, the theoretical framework and conclusions for the positive
mass theorem are examined, including the origin of the total energy
definition, which is usually omitted in the original research works for
various proofs of the theorem r10 ; r11 ; r12 ; r13 ; r14 .

In general relativity, the metric $g_{\mu\nu}$ are the _solutions_ of Einstein
field equation, while the energy–momentum tensor $T^{\mu\nu}$, which causes
space to curve (r1, , p. 5), is the _source_ of the field equation; just like
the EM fields are the solutions of Maxwell equations, while the charge and
current are the source of Maxwell equations (r4, , p. 238). (Note that
Arnowitt-Deser-Misner practice of the tensor indices $\mu,\nu=0,1,2,3$ is used
in this section for the convenience of readers to check with the source papers
r17 ; r18 .)

According to Arnowitt, Deser and Misner r17 , the definition of the total
energy–momentum four-vector $P^{\mu}$ is given by the _volume_ integral of the
components of $T^{0\mu}$, namely

$\displaystyle P^{\mu}=\int T^{0\mu}\textrm{d}^{3}x,$ (1)

which can be expressed as surface integral through Einstein field equations
and Gauss theorem (r1, , p. 462). $P^{0}=E$ in above Eq. (1) is the total
energy r18 .

In the proofs of the positive mass theorem r11 ; r12 ; r14 , the total energy
follows the definition in r17 ; r18  (r1, , p. 462), given by

$\displaystyle E=\int
T^{00}\textrm{d}^{3}x~{}~{}=\frac{1}{16\pi}\int(g_{jk,k}-g_{kk,j})\mathrm{d}^{2}S^{j},$
(2)

where the volume integral is evaluated over the source, and the surface
integral is done over a closed surface completely surrounding the source in
the asymptotically flat region. $T^{00}$ satisfies the dominant energy
condition $T^{00}\geq|T^{\mu\nu}|$ r9 ; r11 ; r14  and it is of the source
term of Einstein field equation, while $(g_{jk,k}-g_{kk,j})$ is the solution
related term.

It should be repeatedly emphasized that $\int(...)\,\textrm{d}^{3}x$ in Eq.
(1) means to carry out _volume_ integration so that Gauss theorem can be used
to convert $\int(...)\,\textrm{d}^{3}x$ into $\int(...)\,\mathrm{d}^{2}S^{j}$
in Eq. (2); confer Ref. (r1, , p. 462).

Arnowitt, Deser and Misner claim that because of the conservation law
$T^{\mu\nu}_{~{}~{}~{},\,\mu}=0$, “$P^{\mu}$ should transform as a four-
vector” r18 , which is clearly endorsed by Nester r19 . In the textbook by
Misner, Thorne and Wheeler (r1, , p. 443, p. 462), it is also emphasized that
the conservation law $T^{\mu\nu}_{~{}~{}\,,\nu}=0$
($T^{\mu\nu}_{~{}~{}~{},\,\mu}=0$) makes $P^{\mu}=\int T^{\mu
0}\textrm{d}^{3}x$ $~{}(P^{\mu}=\int T^{0\mu}\textrm{d}^{3}x)$ be a constant
four-vector.

Thus the positive mass theorem actually states that if the tensor $T^{\mu\nu}$
satisfies the conservation law and dominant energy condition, namely
$T^{\mu\nu}_{~{}~{}~{},\,\mu}=0$ and $T^{00}\geq|T^{\mu\nu}|$ hold, then two
conclusions can be drawn: Conclusion (a) $P^{\mu}=\int
T^{0\mu}\textrm{d}^{3}x$ is a constant four-vector, and Conclusion (b)
$E=P^{0}=\int T^{00}\textrm{d}^{3}x$ is nonnegative.

It is worthwhile to point out that in Eq. (2), no matter in what
asymptotically flat region the solution-related surface integral
$(1/16\pi)\int(g_{jk,k}-g_{kk,j})\mathrm{d}^{2}S^{j}$ is evaluated, it is
always equal to the source-related volume integral $\int
T^{00}\textrm{d}^{3}x$; in other words, asymptotically flat solutions
$g_{\mu\nu}$ of Einstein field equation always satisfies Eq. (2). Thus if
$\int T^{00}\textrm{d}^{3}x$ is not a component of four-vector, then
$(1/16\pi)\int(g_{jk,k}-g_{kk,j})\mathrm{d}^{2}S^{j}$ is not either.

In the following section, we will prove by enumerating a counterexample that
$E=P^{0}=\int T^{00}\textrm{d}^{3}x$ is not a component of four-vector, and
thus the above Conclusion (a) in the positive mass theorem is not true in
general.

##  III Counterexample

In this section, a counterexample of the positive mass theorem is provided,
which is constructed from the wave four-vector of a plane wave of light in
free space. This wave four-vector, first formulated by Einstein in 1905,
predicts the relativistic Doppler effect r16 , which is the physical basis of
the two successive frequency upshifts of a free-electron laser r20  and has
been widely demonstrated by experiments r21 ; r22 ; r23 ; r24 .

Suppose that observed in an inertial frame $XYZ$ in free space, the wave four-
vector for a plane light wave is given by r16

$\displaystyle K^{\mu}=\left(\mathbf{k}_{\textrm{w}},\frac{\omega}{c}\right),$
(3)

where $K^{1,2,3}=(\mathbf{k}_{\textrm{w}})_{x,y,z}$, $K^{4}=\omega/c$,
$\mathbf{k}_{\textrm{w}}$ is the wave vector, $\omega~{}(>0)$ is the angular
frequency, and $c$ is the speed of light in free space.

The wave four-vector satisfies $K^{\mu}K_{\mu}=0$, and we have
$(\omega/c)^{2}-\mathbf{k}_{\textrm{w}}^{2}=0\Rightarrow|\mathbf{k}_{\textrm{w}}|=\omega/c=K^{4}$.
With
$K^{4}=|\mathbf{k}_{\textrm{w}}|=[(\mathbf{k}_{\textrm{w}})_{x}^{2}+(\mathbf{k}_{\textrm{w}})_{y}^{2}+(\mathbf{k}_{\textrm{w}})_{z}^{2}]^{1/2}$
and $(\mathbf{k}_{\textrm{w}})_{x,y,z}=K^{1,2,3}$ taken into account, we have

$K^{4}\geq|K^{1,2,3}|.$ (4)

As a counterexample of the positive mass theorem, the wave four-tensor is
defined as

$T^{\mu\nu}=K^{\mu}K^{\nu}~{}\textrm{for}~{}\mathbf{x}\in
V;\quad=0~{}\textrm{for}~{}\mathbf{x}\not\in V;$ (5)

where the finite volume $V$ $(\neq 0)$ is fixed in $XYZ$. According to above
Eq. (5), the wave four-tensor $T^{\mu\nu}$ is a finite distribution source of
Einstein field equation. Because $K^{\mu}$ is independent of space and time
variables $X^{\mu}=(\mathbf{x},ct)$, $T^{\mu\nu}$ is divergence-free within
the volume $V$, namely

$T^{\mu\nu}_{~{}~{}\,,\nu}=\frac{\partial T^{\mu\nu}}{\partial
X^{\nu}}=\frac{\partial(K^{\mu}K^{\nu})}{\partial X^{\nu}}=0.\\\ $ (6)

From Eqs. (4) and (5) we find that $T^{\mu\nu}$ satisfies the dominant energy
condition in $V$,

$T^{44}=(K^{4})^{2}\geq|K^{\mu}K^{\nu}|=|T^{\mu\nu}|.$ (7)

Consider the volume integral of the components of $T^{\mu 4}$ over $V$, given
by

$\displaystyle P^{\mu}$ $\displaystyle=\int_{V}T^{\mu 4}\textrm{d}^{3}x$
$\displaystyle=\int_{V}(K^{\mu}K^{4})\textrm{d}^{3}x=(K^{4}V)K^{\mu},$ (8)

where $\int_{V}(K^{\mu}K^{4})\,\textrm{d}^{3}x=(K^{4}V)K^{\mu}$ is employed
because $K^{\mu}$ is independent of $X^{\mu}=(\mathbf{x},ct)$.

From above Eqs. (6) and (7), we know that the wave four-tensor $T^{\mu\nu}$
satisfies the conservation law $T^{\mu\nu}_{~{}~{}\,,\nu}=0$ and the dominant
energy condition $T^{44}\geq|T^{\mu\nu}|$ over $V$. Thus according to the
positive mass theorem, the quantity $P^{\mu}=(K^{4}V)K^{\mu}$, given by Eq.
(8), must be a Lorentz covariant four-vector.

Thus if we can prove that $P^{\mu}=(K^{4}V)K^{\mu}$ is not a Lorentz four-
vector, then $T^{\mu\nu}$ becomes a counterexample of the positive mass
theorem. As shown below, $P^{\mu}=(K^{4}V)K^{\mu}$ indeed cannot be a four-
vector.

Suppose that the inertial frame $X^{\prime}Y^{\prime}Z^{\prime}$ moves with
respect to $XYZ$ at an arbitrary velocity of $\mathbf{v}=\boldsymbol{\beta}c$,
with $X^{\prime}Y^{\prime}Z^{\prime}$ and $XYZ$ overlapping at
$t^{\prime}=t=0$ r6 . Observed in $X^{\prime}Y^{\prime}Z^{\prime}$, the volume
integral of the components of $T^{\prime\mu 4}$ is given by

$\displaystyle P^{\prime\mu}$ $\displaystyle=\int_{V^{\prime}}T^{\prime\mu
4}\textrm{d}^{3}x^{\prime}$
$\displaystyle=\int_{V^{\prime}}(K^{\prime\mu}K^{\prime
4})\textrm{d}^{3}x^{\prime}=(K^{\prime 4}V^{\prime})K^{\prime\mu},$ (9)

where the volume $V^{\prime}$ is moving with respect to
$X^{\prime}Y^{\prime}Z^{\prime}$ at the velocity
$\mathbf{v}^{\prime}=\boldsymbol{\beta}^{\prime}c=-\boldsymbol{\beta}c$.

Comparing above Eq. (9) with Eq. (8), we find that if the relation between
$P^{\mu}=(K^{4}V)K^{\mu}$ and $P^{\prime\mu}=(K^{\prime
4}V^{\prime})K^{\prime\mu}$ follows four-vector Lorentz transformation, then
$(K^{4}V)=(K^{\prime 4}V^{\prime})$ = _Lorentz covariant scalar_ must hold,
because $K^{\mu}$ is a four-vector. However as shown below,
$(K^{4}V)=(K^{\prime 4}V^{\prime})$ cannot hold for an _arbitrary_
$\mathbf{v}=\boldsymbol{\beta}c$, namely $(K^{4}V)$ is not a covariant scalar,
and so $P^{\mu}=(K^{4}V)K^{\mu}$ cannot be a four-vector.

According to the _change of variables theorem_ in principles of classical
mathematical analysis (r25, , p. 252), the transformation of differential
elements appearing in Eq. (8) and Eq. (9) is given by

$\textrm{d}^{3}x^{\prime}=\left|\frac{\partial(x^{\prime},y^{\prime}z^{\prime})}{\partial(x,y,z)}\right|\textrm{d}^{3}x=\frac{1}{\gamma}\textrm{d}^{3}x,$
(10)

where $\gamma=(1-\boldsymbol{\beta}^{2})^{-1/2}$ is the Lorentz factor, and
the Jacobi determinant
$\partial(x^{\prime},y^{\prime}z^{\prime})/\partial(x,y,z)=1/\gamma$ is
employed. It should be noted that the transformation of above Eq. (10) meets
the physical requirement of Lorentz contraction effect of the length of a
moving rigid rod argued by Einstein r16  and further formulated in Ref. (r26,
, Footnote 2 there).

From above Eq. (10), we have

$\int_{V^{\prime}}\textrm{d}^{3}x^{\prime}=\int_{V}\left(\frac{1}{\gamma}\right)\textrm{d}^{3}x\quad\quad\textrm{or}\quad\quad
V^{\prime}=\frac{V}{\gamma}.$ (11)

According to Einstein r16 , the Lorentz transformation between $K^{\prime
4}=\omega^{\prime}/c$ and $K^{4}=\omega/c$ is given by

$K^{\prime 4}=\gamma(K^{4}-\mathbf{k}_{\textrm{w}}\cdot\boldsymbol{\beta}).$
(12)

Combining Eqs. (11) and (12), we have

$K^{\prime
4}V^{\prime}=K^{4}V-(\mathbf{k}_{\textrm{w}}\cdot\boldsymbol{\beta})V,$ (13)

which is valid for any $|\boldsymbol{\beta}|<1$.

Obviously, $(\mathbf{k}_{\textrm{w}}\cdot\boldsymbol{\beta})V=0$ cannot hold
for arbitrary $|\boldsymbol{\beta}|<1$; for example,
$(\mathbf{k}_{\textrm{w}}\cdot\boldsymbol{\beta})V=(\omega/c)|\boldsymbol{\beta}|V\neq
0$ holds for $\mathbf{k}_{\textrm{w}}\|\,\boldsymbol{\beta}$ and
$\boldsymbol{\beta}\neq 0$. Thus from Eq. (13) we conclude that
$(K^{4}V)=(K^{\prime 4}V^{\prime})$ cannot hold for arbitrary
$|\boldsymbol{\beta}|<1$, namely $(K^{4}V)$ is not a Lorentz covariant scalar.
Since $(K^{4}V)$ is not a covariant scalar, the quantity
$P^{\mu}=(K^{4}V)K^{\mu}$, given by Eq. (8), cannot be a four-vector.

So far, we have finished the proof that the wave four-tensor $T^{\mu\nu}$
satisfies $T^{\mu\nu}_{~{}~{}\,,\nu}=0$ and $T^{44}\geq|T^{\mu\nu}|$ over the
volume $V$, but $P^{\mu}=\int_{V}T^{\mu 4}\textrm{d}^{3}x$ is not a four-
vector ($\Leftrightarrow$ $P^{4}=\int_{V}T^{44}\textrm{d}^{3}x$ is not a
component of four-vector); namely the wave four-tensor $T^{\mu\nu}$, defined
in Eq. (5), is a counterexample of the positive mass theorem.

##  IV Conclusions and remarks

By examining the background, we have found that the positive mass theorem has:

  * •

_two requirements_ : the momentum–energy tensor $T^{\mu\nu}$, as the source of
Einstein field equation, satisfies (a) the conservation law
$T^{\mu\nu}_{~{}~{}\,,\nu}=0$ and (b) the dominant energy condition
$T^{44}\geq|T^{\mu\nu}|$;

  * •

_two conclusions_ : (a) the time-column volume integral $P^{\mu}=\int T^{\mu
4}\textrm{d}^{3}x$ constitutes a total four-momentum that is conserved (r1, ,
p. 443), and (b) the total energy $E=P^{4}=\int T^{44}\textrm{d}^{3}x$ is
nonnegative.

We have proved that the wave four-tensor, given in Eq. (5), is a
counterexample of the positive mass theorem because it satisfies the above two
requirements, but Conclusion (a) does not hold. Since Conclusion (a) is not
valid, Conclusion (b) naturally loses its physical significance. Thus the
positive mass theorem is fundamentally flawed.

From the counterexample provided in the present paper, one can see that the
problem of positive mass theorem comes from the traditional conservation law
$T^{\mu\nu}_{~{}~{}\,,\nu}=0$, which is actually an incorrect conjecture r5 ;
r6 ; r7 ; r8  although some proofs in support of it are provided in textbooks,
such as in (r1, , p. 142) (r27, , p. 318).

In addition to the conservation law $T^{\mu\nu}_{~{}~{}\,,\nu}=0$, there are
some other basic concepts and definitions which have not yet reached a
consensus, such as the covariance of physical laws r28 , electromagnetic power
flow in an anisotropic medium r29 , and Fermat principle for a plane light
wave r30 . Probably one of the best known and most controversial examples is
about how to define a particle’s mass. In respected textbooks r15 ; r31 , two
different mass definitions are argued, which leads to an amazing and
unexpected result: the Lorentz covariance of Dirac wave equation is not
compatible with Einstein mass-energy equivalence, as shown below.

In the first mass definition for a particle, the four-momentum is given by
(r15, , p. 110)

$\displaystyle P^{\mu}=\frac{m}{\gamma}U^{\mu}$
$\displaystyle=\left(m\mathbf{v},~{}\frac{mc^{2}}{c}\right)$
$\displaystyle=\left(\mathbf{p},~{}\frac{E}{c}\right),\quad(\textrm{First
definition})$ (14)

where $m$ is defined as the particle’s (relativistic inertial) mass;
$(m/\gamma)$ is a _covariant_ scalar (invariant) r28 ;
$U^{\mu}=\gamma(\mathbf{v},c)$ is the four-velocity, with $\mathbf{v}$ the
particle’s velocity; $\mathbf{p}$ is the momentum; and $E$ is the energy. In
this definition, (i) $(m/\gamma)=m_{0}$ holds, where $m_{0}$ is the rest mass
and it is $not$ a covariant invariant r28 ; (ii) Einstein mass–energy
equivalence equation $E=mc^{2}$ holds in all inertial frames.

In the second mass definition, the particle’s four-momentum is given by (r31,
, p. 289)

$\displaystyle P^{\mu}=mU^{\mu}$ $\displaystyle=\left(\gamma
m\mathbf{v},~{}\frac{\gamma mc^{2}}{c}\right)$
$\displaystyle=\left(\mathbf{p},~{}\frac{E}{c}\right),\quad(\textrm{Second
definition})$ (15)

where $m$ is defined as the particle’s (invariant) mass and it is a covariant
invariant; $\mathbf{p}=\gamma m\mathbf{v}$ is the momentum; and
$E=\gamma(mc^{2})$ is the energy. In this definition, (i) $m=m_{0}$ always
holds in all inertial frames (namely $m_{0}$ and $m$ are the same), and
$m_{0}$ is a covariant invariant; (ii) $E=\gamma(mc^{2})$ holds, but the
Einstein mass–energy equivalence equation $E=(mc^{2})$ does not hold except
for in the particle-rest frame.

According to the first mass definition shown in Eq. (14), the particle’s rest
mass is not a covariant invariant; thus Dirac wave equation for an electron is
not covariant because the electron’s rest mass appearing in Dirac equation is
not a covariant invariant although Dirac took it to be in his proof (r32, , p.
258). However the covariance of Dirac equation is consistent with the second
definition shown in Eq. (15), which contradicts with Einstein mass–energy
equivalence equation. From this one can conclude that the Lorentz covariance
of Dirac wave equation is not compatible with Einstein mass–energy
equivalence. In other words, if Einstein mass–energy equivalence is valid,
then Dirac wave equation is not Lorentz covariant, and vice versa.

It is worthwhile to point out that in the first mass definition, the photon
has its (kinetic) mass r28  because “the mass of a body is a measure of its
energy content” r33 , while in the second mass definition, the photon does not
have any mass (equal to zero) (r34, , p. 99). However Einstein’s thought
experiment for mass–energy equivalence does support the view that the photon
has a non-zero mass, as shown below.

In his thought experiment r33 , Einstein assumes that a stationary body with a
rest energy of $E_{0}$ emits two identical photons (“plane waves of light” in
Einstein’s words) in free space at the same time in the opposite directions so
that the body keeps at rest after the emissions. According to the energy
conservation law, Einstein argues

$E_{0}=(E_{0}-\Delta E_{0})+\frac{1}{2}\Delta E_{0}+\frac{1}{2}\Delta E_{0},$
(16)

where $(E_{0}-\Delta E_{0})$ is the body’s energy after emissions, and the
last two terms of $(\Delta E_{0}/2)$ are the (kinetic) energies of the two
photons respectively. Dividing above Eq. (16) by $c^{2}$ yields

$\frac{E_{0}}{c^{2}}=\left(\frac{E_{0}}{c^{2}}-\frac{\Delta
E_{0}}{c^{2}}\right)+\frac{1}{2}\frac{\Delta
E_{0}}{c^{2}}+\frac{1}{2}\frac{\Delta E_{0}}{c^{2}}.$ (17)

In terms of the principle of relativity, Einstein proved that the body’s rest
mass is reduced by $(\Delta E_{0}/c^{2})$ after emissions. Thus it is
legitimate to define $(E_{0}/c^{2})$ and $(E_{0}-\Delta E_{0})/c^{2}$ as the
masses of the stationary body before and after the emissions. From this it
follows that the last two terms of $\Delta E_{0}/(2c^{2})$ in Eq. (17)
legitimately denote the masses of the two photons respectively.

Thus Einstein’s thought experiment r33  supports the conclusion that the
photon has a non-zero kinetic energy, and it has a non-zero mass; namely the
mass and energy are equivalent, as Einstein argued.

## References

  * (1) C.W. Misner, K.S. Thorne, and A.Wheeler, Gravitation, W.H. Freeman and Company, San Francisco, 1973.
  * (2) S. Weinberg, Gravitation and Cosmology: Principles and Applications of the General Theory of Relativity, John Wiley & Sons, New York, 1972.
  * (3) C. Møller, The Theory of Relativity, Oxford University Press, London, 1955.
  * (4) J.D. Jackson, Classical Electrodynamics, 3rd Edition, John Wiley & Sons, NJ, 1999.
  * (5) C. Wang, von Laue’s theorem and its applications, Can. J. Phys. 93 (2015) 1470–1476.
  * (6) C. Wang, Self-consistent theory for a plane wave in a moving medium and light-momentum criterion, Can. J. Phys. 93 (2015) 1510–1522.
  * (7) C. Wang, Disproof of a widely-accepted mathematical conjecture, Optik 140 (2017) 1110–1113; arXiv:1612.03201.
  * (8) C. Wang, Resolution of two fundamental issues in the dynamics of relativity and exposure of a real version of the emperor’s new clothes, arXiv:1702.03200.
  * (9) S. Hawking, The conservation of matter in general relativity, Commun. Math. Phys. 18 (1970) 301–306.
  * (10) R. Schoen and S.T. Yau, On the proof of the positive mass conjecture in general relativity, Commun. Math. Phys. 65 (1979) 45–76.
  * (11) R. Schoen and S.T. Yau, Proof of the positive mass theorem. II, Commun. Math. Phys. 79 (1981) 231–260.
  * (12) E. Witten, A new proof of the positive energy theorem, Commun. Math. Phys. 80 (1981) 381–402.
  * (13) R. Schoen and S.T. Yau, Positivity of the total mass of a general space-time, Phys. Rev. Lett. 43 (1979) 1457-1459.
  * (14) T. Parker and C.H. Taubes, On Witten’s proof of the positive energy theorem, Commun. Math. Phys. 84 (1982) 223–238.
  * (15) W. Rindler, Relativity: Special, General, and Cosmological, 2nd Edition, Oxford, NY, 2006.
  * (16) A. Einstein, Zur Elektrodynamik bewegter Körper, Ann. Phys. 322 (1905) 891–921.
  * (17) R. Arnowitt, S. Deser, and C.W. Misner, Energy and the criteria for radiation in general relativity, Phys. Rev. 118 (1960) 1100–1104.
  * (18) R. Arnowitt, S. Deser, and C. W. Misner, Coordinate invariance and energy expressions in general relativity, Phys. Rev. 122 (1961) 997–1006.
  * (19) J.M. Nester, A new gravitational energy expression with a simple positivity proof, Phys. Lett. A 83 (1981) 241–242.
  * (20) C. Wang, The relativistic Doppler effect: when a zero-frequency shift or a red shift exists for sources approaching the observer, Ann. Phys. (Berlin) 523 (2011) 239–246. As shown in this reference, after two successive frequency upshifts in a magnetic wiggler free-electron laser (FEL), the radiation wavelength is given by $\lambda=\lambda_{\mathrm{w}}(1-\beta_{\|})/\beta_{\|}$, where $\lambda_{\mathrm{w}}$ is the wiggler period, and $\beta_{\|}$ is the electron’s axial drift velocity normalized to the vacuum light speed $c$. When a high-energy and small wiggler-velocity electron beam is employed so that $(1+a_{\textrm{w}}^{2})/\gamma^{2}<<1$ holds, the radiation wavelength can be well approximated as $\lambda=\lambda_{\mathrm{w}}(1+a^{2}_{\mathrm{w}})/(2\gamma^{2})$ for a helical wiggler FEL r21 ; r22 , and $\lambda=\lambda_{\mathrm{w}}(1+a^{2}_{\mathrm{w}}/2)/(2\gamma^{2})$ for a planar wiggler FEL r23 ; r24 , where $\gamma$ is the electron’s relativistic factor; $a_{\mathrm{w}}=(|e|B_{\mathrm{w}}/m)(\lambda_{\mathrm{w}}/2\pi c)=0.934\,B_{\mathrm{w}}[\mathrm{T}]\,\lambda_{\mathrm{w}}[\mathrm{cm}]$ is the wiggler parameter, with $e$ and $m$, respectively, the electron’s charge and rest mass, and $B_{\mathrm{w}}$ the amplitude of wiggler magnetic field. (1) For a helical wiggler FEL with $\gamma=47.97$ (24 MeV), $\lambda_{\mathrm{w}}=3.2$ cm, and $B_{\mathrm{w}}=0.24$ T, we have $a_{\mathrm{w}}=0.717$ and the theoretical radiation wavelength is given by $\lambda=3.2\,\mathrm{cm}\times(1+0.717^{2})/(2\times 47.97^{2})=10.53$ $\upmu$m, compared with the experimental observation of 10.6 $\upmu$m in Ref. r21 . (2) For a planar wiggler FEL with $\gamma=130.2$ (66 MeV), $\lambda_{\mathrm{w}}=3.6$ cm, and $B_{\mathrm{w}}=0.29$ T, we have $a_{\mathrm{w}}=0.975$ and the theoretical radiation wavelength is given by $\lambda=3.6\,\mathrm{cm}\times(1+0.975^{2}/2)/(2\times 130.2^{2})=1.567$ $\upmu$m, compared with the experimental observation of 1.57 $\upmu$m in Ref. r23 . Thus the free-electron laser has provided strong and convincing experimental support for the wave four-vector, which is one of the basic results of Einstein’s special relativity r16 .
  * (21) L.R. Elias, W.M. Fairbank, J.M.J. Madey, H.A. Schwettman, and T.I. Smith, Observation of stimulated emission of radiation by relativistic electrons in a spatially periodic transverse magnetic field, Phys. Rev. Lett. 36 (1976) 717–720.
  * (22) D.A.G. Deacon, L.R. Elias, J.M.J. Madey, G.J. Ramian, H.A. Schwettman, and T.I. Smith, First operation of a free-electron laser, Phys. Rev. Lett. 38 (1977) 892–894.
  * (23) J.A. Edighoffer, G.R. Neil, C.E. Hess, T.I. Smith, S.W. Fornaca, and H.A. Schwettman, Variable-wiggler free-electron-laser oscillation, Phys. Rev. Lett. 52 (1984) 344–347.
  * (24) J. Andruszkow, B. Aune, V. Ayvazyan, N. Baboi, R. Bakker, V. Balakin, D. Barni, A. Bazhan, M. Bernard, A. Bosotti _et al._ , First observation of self-amplified spontaneous emission in a free-electron laser at 109 nm wavelength, Phys. Rev. Lett. 85 (2000) 3825–3829.
  * (25) W. Rudin, Principles of Mathematical Analysis, 3rd Edition, McGraw-Hill, New York, 1976.
  * (26) C. Wang, Minkowski tensor in electrodynamics of moving media and three rules for construction of the physical tensors in Einstein’s special relativity, Optik 223 (2020) 165469.
  * (27) W. Thirring, Classical Mathematical Physics, 3rd Edition, Springer-Verlag, New York, 1997.
  * (28) C. Wang, Criterion for testing the covariance of physical laws and Gordon optical metric, Optik 239 (2021) 166775.
  * (29) C. Wang, Electromagnetic power flow, Fermat’s principle, and special theory of relativity, Optik 126 (2015) 2703–2705.
  * (30) C. Wang, New insight into light propagation and light-matter interactions with applications to experimental observations, Optik 204 (2020) 163954.
  * (31) H. Goldstein, C. Poole, and J. Safko, Classical Mechanics, 3rd Edition, Addison-Wesley, NY, 2000.
  * (32) P.A.M. Dirac, The Principles of Quantum Mechanics, 4th Edition, Oxford, London, 1958.
  * (33) A. Einstein, Ist die Trägheit eines Körpers von seinem Energieinhalt abhängig? Ann. Phys. 323 (1905) 639–641.
  * (34) D. Griffiths, Introduction to Elementary Particles, 2nd Edition, Wiley-VCH, Weinheim, 2008.


# Modeling and Control of Diesel Engine Emissions using Multi-layer Neural
Networks and Economic Model Predictive Control

Jiadi Zhang  Xiao Li  Mohammad Reza Amini  Ilya Kolmanovsky  Munechika
Tsutsumi  Hayato Nakada  Department of Aerospace Engineering, University of
Michigan, Ann Arbor, MI 48109, USA (e-mails<EMAIL_ADDRESS>Department of Naval Architecture and Marine Engineering, University of
Michigan, Ann Arbor, MI 48109, USA (e-mail<EMAIL_ADDRESS> Hino Motors,
Ltd., Tokyo 191-8660, Japan (e-mails<EMAIL_ADDRESS>
###### Abstract

This paper presents the results of developing a multi-layer Neural Network
(NN) to represent diesel engine emissions and integrating this NN into control
design. Firstly, a NN is trained and validated to simultaneously predict
oxides of nitrogen ($NOx$) and $Soot$ using both transient and steady-state
data. Based on the input-output correlation analysis, inputs to NN with the
highest influence on the emissions are selected while keeping the NN structure
simple. Secondly, a co-simulation framework is implemented to integrate the NN
emissions model with a model of a diesel engine airpath system built in GT-
Power and used to identify a low-order linear parameter-varying (LPV) model
for emissions prediction. Finally, an economic supervisory model predictive
controller (MPC) is developed using the LPV emissions model to adjust
setpoints to an inner-loop airpath tracking MPC. Simulation results are
reported illustrating the capability of the resulting controller to reduce
$NOx$, meet the target $Soot$ limit, and track the adjusted intake manifold
pressure and exhaust gas recirculation (EGR) rate targets.

###### keywords:

Model predictive control, Linear parameter-varying control, Neural networks,
Diesel engine control

, , , , ,

##  1 Introduction

To achieve emission reduction in diesel engines, advanced engine control
methods can be leveraged to coordinate engine fueling, exhaust gas
recirculation (EGR) valve, and variable geometry turbocharger (VGT). In
particular, Model Predictive Control (MPC) has shown significant promise for
improving engine transient control while satisfying state and control
constraints (Norouzi et al. (2021)). Accurate and low-complexity engine-out
emission models can support control design and verification before physical
engine testing (Li et al. (2017)).

Different emission modeling methods have been proposed in the literature
including map-based quasi-steady models (Hagena et al. (2006)), physics-based
grey-box models (Heywood (2018); Arrègle et al. (2008)), and machine learning
based models (Winkler-Ebner et al. (2010); Prokhorov (2008); De Cesare and
Covassin (2011); Li et al. (2017)). As shown in Hagena et al. (2006), map-
based models may not accurately predict transient effects such as smoke spikes
during fuel tip-ins. On the other hand, physics-based grey-box models often
exploit in-cylinder pressure as an input, which is not available in commercial
vehicles in production. Hence, in this paper, a data-driven machine learning-
based model is developed using limited training data to support the emissions-
oriented MPC design.

The existing literature on the application of MPC to diesel engines has
addressed intake manifold pressure and airflow/EGR rate setpoint tracking in
the airpath system (Huang et al. (2013, 2018); Ortner and Del Re (2007)) as
well as higher-level objectives such as enhanced fuel economy or reduced
emissions (Huang et al. (2020)). Liao-McPherson et al. (2020); Liu et al.
(2021) have developed solutions along the lines of nonlinear economic MPC
(EMPC) to control $NOx$ and $Soot$ emissions. Broomhead et al. (2016) have
used EMPC for power tracking while limiting emissions in diesel generator
applications.

In this paper we focus on a hierarchical MPC architecture for emissions
reduction with an outer-loop (supervisory) economic MPC which adjusts fueling
rate and intake manifold pressure and EGR rate targets (set-points) to an
inner-loop tracking MPC for the airpath system. For the inner-loop MPC, we
adopt the airpath controller design from Liao-McPherson et al. (2020); Zhang
et al. (2022a), which exploits feedforward MPC and feedback rate-based MPC.
Note that in automotive applications, a feedback tracking controller is often
coupled with a feedforward controller to speed up the transient response
(Norouzi et al. (2021)). One of the most common feedforward controllers is
based on look-up tables for set-points and actuator positions, however, an
MPC-based feedforward (Liao-McPherson et al. (2020); Zhang et al. (2022a)) has
been shown to improve transient response.

To support controller implementation, we first develop a multi-layered Neural
Network (NN) model for engine-out $NOx$ and $Soot$ emissions based on
experimental dynamometer data. Then, we integrate this NN-based emissions
model with a high-fidelity engine model in GT-Power that predicts flows,
pressures, and temperatures in the engine to obtain a high-fidelity model for
closed-loop simulations. We then use this high-fidelity model to identify a
low-order linear parameter-varying (LPV) model for emissions prediction, and
we use this LPV model as a prediction model to implement supervisory economic
MPC. Closed-loop simulation results are finally reported followed by
concluding remarks.

##  2 Diesel Engine Modeling

The schematic of the diesel engine considered in this paper is shown in Figure
1. The EGR system dilutes the cylinder charge and lowers the peak combustion
temperature, thereby reducing oxides of nitrogen $NOx$ emissions. The flow
from the exhaust manifold of the engine to the intake manifold is controlled
by the EGR valve. As the EGR flow also depends on the exhaust pressure, this
flow is also affected by the VGT. Both VGT and EGR valve affect the intake
manifold pressure and air-to-fuel ratio. High levels of recirculated exhaust
gas and low air-to-fuel engine can lead to increased smoke/soot emissions. The
engine operating point is typically defined by the engine speed ($N_{\tt e}$)
and total fuel injection rate ($w_{\tt inj}$) inferred from the engine torque
request; these are treated as the exogenous inputs.

Figure 1: Diesel engine with exhaust gas recirculation (EGR) and a variable
geometry turbocharger (VGT).

A typical diesel engine airpath controller tracks target values (set-points)
for the intake manifold pressure ($p_{\tt im}$) and EGR rate ($\chi_{\tt
egr}$). The EGR rate is defined by

$\chi_{\tt egr}=\frac{w_{\tt egr}}{w_{\tt egr}+w_{\tt c}},$ (1)

where $w_{\tt c}$ is the mass flow into the intake manifold through the
compressor and intercooler, and $w_{\tt egr}$ is the mass flow through the EGR
valve from the exhaust manifold into the intake manifold.

In this paper, the objective is to coordinate the EGR valve and VGT actuators
to control the intake manifold pressure ($p_{\tt im}$) and EGR rate
($\chi_{\tt egr}$) to target values that are computed and provided to the
airpath controller by a supervisory controller to satisfy emissions and fuel
consumption requirements.

As in our previous work on diesel engine airpath control (Zhang et al. (2022b,
a)), our controller development relies on a high-fidelity diesel engine model
in GT-Power. This model represents the responses of flows, and pressures in
different parts of the engine at a crankangle resolution and has been
validated against experimental data from engine dynamometer testing. This
model, however, does not represent engine feedgas emissions; hence, to enable
emissions-oriented controller development in this paper, we first augment this
GT-Power engine model with data-driven emissions models. We then develop
control-oriented models for feedgas emissions and airpath, and use them for
EMPC design.

###  2.1 Data-driven Engine Feedgas Emissions Modeling

Neural Networks, as universal function approximators, are capable of learning
emission models from data. Given experimental dynamometer data sets with both
transient and steady-state measurements, we model the diesel engine emissions
using a multi-layered NN (Figure 2), where the output of each layer is defined
by

$y_{i}=\sigma_{ReLU}\left(W_{i}y_{i-1}+b_{i}\right),~{}~{}i=1,\dots,4,$ (2)

and where $W_{i}$ and $b_{i}$ are the network parameters of the $i$th layer
while $y_{0}\in\mathbb{R}^{10}$ and $y_{4}\in\mathbb{R}^{2}$ are the inputs
and predicted emissions outputs, respectively. The NN is trained on an
experimental dataset consisting of both steady-state ($\tt ss$)
$\\{y_{0,i}^{\tt ss},y_{4,i}^{\tt ss}\\}_{i=1,\dots,306}$ and transient ($\tt
ts$) measurements $\\{y_{0,i}^{\tt ts},y_{4,i}^{\tt ts}\\}_{i=1,\dots,12001}$.
The objective of the training a NN model is to minimize the sums of squares of
the errors between the output $y_{4,i}$, $i=1,2$, and the actual emissions,
namely, oxides of nitrogen ($NOx$) and $Soot$.

Figure 2: Schematic of the developed multi-layer neural network diesel
emissions model.

The three hidden layers of our NN model of size $1024$, $512$, and $32$,
respectively, are chosen consistently with an encoder-decoder architecture
while keeping the number of trainable parameters relatively small. The first
hidden layer, functioning as an encoder, blends the input measurements from
$\mathbb{R}^{10}$ into features in a higher dimension $\mathbb{R}^{1024}$. The
consecutive two hidden layers have $512$ and $32$ neurons, respectively. They
decode the high-dimensional features and project the vector layer-by-layer
into $\mathbb{R}^{2}$.

Table 1: Normalized cross-covariance with zero lag between measured input variable candidates and emissions. | Transient | Steady State
---|---|---
Measured | $NOx$ | $Soot$ | $NOx$ | $Soot$
Variables | $\rm[ppm]$ | $\rm[\%]$ | $\rm[ppm]$ | $\rm[\%]$
Injection Pressure $\rm[MPa]$ | 0.45 | -0.32 | -0.09 | -0.18
Main injection | 0.27 | -0.39 | -0.03 | -0.19
timing $\rm[BTDC]$
Main injection fuel | 0.68 | -0.12 | 0.38 | -3.0e-3
flow rate $\rm[mm^{3}/st]$
Pre-injection fuel | 0.22 | 0.07 | -8.2e-17 | -1.7e-16
flow rate $\rm[mm^{3}/st]$
Engine torque output $\rm[Nm]$ | 0.67 | -0.10 | 0.38 | -0.02
Engine speed $\rm[rpm]$ | 0.28 | -0.40 | -0.25 | -0.07
Intake manifold | 0.47 | -0.21 | 0.15 | 0.15
pressure $\rm[kPa]$
Exhaust manifold | 0.49 | -0.22 | 0.11 | -0.04
pressure $\rm[kPa]$
Mass air flow $\rm[G/s]$ | 0.39 | -0.27 | 0.04 | -0.12
EGR position $\rm[\%]$ | -0.23 | -0.33 | 0.10 | 0.22
VGT position $\rm[\%]$ | -0.32 | 0.42 | 0.08 | 0.13

To enable our NN model to exploit two datasets (i.e., steady-state dataset and
transient dataset) of distinct physical nature and dataset sizes we used the
following procedure. Firstly, the input variables used in training were
selected based on correlation analysis. For instance, Table 1 indicates low
normalized cross-covariance for the pre-injection fuel flow rate; hence we
remove it from the set of model inputs while other signals in Table 1 are kept
in $y_{0}\in\mathbb{R}^{10}$. Secondly, we exclude outliers from the transient
training data set to facilitate the integration of the transient and steady
state data. For this, we empirically choose a threshold, $\epsilon>0$, and
exclude a transient data point $(y_{0,i}^{\tt ts},y_{4,i}^{\tt ts})$ if

$r(y_{0,i}^{\tt ts})=\left((y_{0,i}^{\tt
ts}-\mu_{ss})^{T}\Sigma_{ss}^{-1}(y_{0,i}^{\tt
ts}-\mu_{ss})\right)^{1/2}>\epsilon,$ (3)

where $r(y_{0,i}^{\tt ts})$ designates the Mahalanobis distance, while
$\mu_{ss}$ and $\Sigma_{ss}$ stand for the mean and the covariance matrix of
$\\{y_{0,i}^{\tt ss}\\}_{i=1,\dots,306}$. Thirdly, note that the transient
dataset of $12001$ data points is significantly larger than the steady-state
dataset of $306$ data points. To avoid overemphasizing transient response in
training the NN model, the final training dataset comprises the transient
dataset after outlier removal and the steady state dataset being duplicated
sevenfold.

Using the Pytorch package developed by Paszke et al. (2019), the NN model was
trained for $1000$ epochs with MSE-loss function and a batch size of $40$. The
stochastic gradient descent algorithm with momentum $\rho=0.9$ was used for
training. The initial learning rate is set to $lr=10^{-4}$ with a decay rate
of $\gamma=0.5$ for every $100$ epochs. The number of data points used for
training, validation, and testing were chosen in the ratio of 70%:15%:15%.

Figure 3 compares the prediction results using the NN model versus the actual
measured emissions data from the testing fraction of the dynamometer dataset.
In Figure 3 and subsequent figures, we are not able to report the $y$-axis
values in order to protect OEM proprietary data. Compared to the measured
$NOx$ data, the predicted $NOx$ has mean errors of 49.8 $\rm ppm$ during
transients and 129.8 $\rm ppm$ in steady state. The mean errors in $Soot$
prediction are 0.45% in transients and 5.7% in steady state. As shown in
Figure 3, our predictions for both transient data (in blue dash lines) and
steady-state data (in blue crosses) match those from the dynamometer data with
relatively small errors on average.

Figure 3: The NN-based diesel emissions prediction performance against the
actual dynamometer data.

###  2.2 Emissions Control-oriented Modeling

To reduce diesel engine emissions, we employ an economic MPC (EMPC). To
implement EMPC, a control-oriented model is needed for predicting the
emissions over the MPC optimization horizon. To that end, in this section, we
develop and validate an LPV model for emissions. To keep the prediction model
as simple as possible, only $NOx$ and $Soot$ are selected as the LPV model
states ($x$). Both emission outputs are assumed to be measured or accurately
estimated. The input ($u$) for our model includes the intake manifold pressure
($p_{\tt im}$), EGR rate ($\chi_{\tt egr}$) and fuel injection rate ($w_{\tt
inj}$).

Both the states and the input variables have been normalized by the steady-
state values ($x^{\tt ss}(\rho)$, $u^{\tt ss}(\rho)$) and their corresponding
standard deviation ($\sigma_{x}(\rho)$, $\sigma_{u}(\rho)$) as follows

$\tilde{x}(\rho)=\frac{x-x^{\tt
ss}(\rho)}{\sigma_{x}(\rho)},\quad\tilde{u}(\rho)=\frac{u-u^{\tt
ss}(\rho)}{\sigma_{u}(\rho)}.$ (4)

The prediction model has an LPV form,

$\tilde{x}_{k+1}(\rho_{k})=\tilde{A}(\rho_{k})\tilde{x}_{k}(\rho_{k})+\tilde{B}(\rho_{k})\tilde{u}_{k}(\rho_{k}),$
(5)

where $k$ denote the discrete time, $\rho_{k}$ is the vector of engine speed
and fuel injection rate at time instant $k$,
$\tilde{A}:\mathbb{R}^{2}\to\mathbb{R}^{2\times 2}$,
$\tilde{B}:\mathbb{R}^{2}\to\mathbb{R}^{2\times 3}$ are operating condition
($\rho$) dependent matrices.

We linearly interpolate the matrices $\tilde{A}(\rho)$ and $\tilde{B}(\rho)$
from $99$ models identified at pre-selected operating points ($\rho$) defined
by 9 values of the engine speed and 11 values of fuel injection rate that
cover the engine operating range. A co-simulation framework is implemented to
integrate the NN emissions model with a model of the diesel engine airpath
system built in GT-Power in order to generate input-output response data
corresponding to small EGR and VGT position perturbations; these data are then
used for local model identification at each of the $99$ operating points. The
data generation procedure is summarized in Figure 4. Using the function ${\tt
n4sid}$ in MATLAB with prediction horizon setting as 50 steps, we could get
the optimized $\tilde{A}$ and $\tilde{B}$ matrix for each operating point that
minimizes the one-step-ahead prediction error between measured and predicted
outputs.

Figure 4: Data generation procedure for identification of control-oriented LPV
model of diesel emissions, where $p_{\tt im}$ is intake manifold pressure,
$\chi_{\tt egr}$ is EGR rate, and $w_{\tt inj}$ is fuel injection quantity.

Figure 5 compares the results of simulating the LPV model (5) and the GT-Power
model with NN-based emissions model over the first $600s$ of the World
Harmonized Transient Cycle (WHTC) driving cycle, where both models respond to
the same input trajectories of engine speed, fuel injection rate, target
intake manifold pressure, and target EGR rate. The trajectories of the target
intake manifold pressure and target EGR rate for these simulations were
generated from look-up tables that specify desired $p_{\tt im}$ and $\chi_{\tt
egr}$ in steady-state as functions of the engine speed and fuel injection
rate.

The average errors in $NOx$ and $Soot$ LPV models as compared to the higher-
fidelity NN model are less than $100~{}ppm$ and 2%, respectively, even though
larger errors are observed in transients. While the transient errors could
potentially be reduced with higher-order LPV or/and NN models, as we will
demonstrate in the following section, the proposed simple LPV model could
successfully support an EMPC design.

Figure 5: The comparison results between control-oriented LPV emissions model
against the GT-Power model co-simulation with NN-based emissions model over
the first 600 $s$ of WHTC driving cycle: (a) $NOx$, and (b) $Soot$.

###  2.3 Diesel Airpath Control-Oriented Model

A control-oriented model of the diesel airpath system is adopted to design the
inner-loop tracking MPC. The airpath model estimates the response of the
intake manifold pressure and EGR rate to EGR and VGT actuators. The intake
manifold pressure ($p_{\tt im}$) and EGR rate ($\chi_{\tt egr}$) are two of
the airpath model states ($z$). The model inputs ($v$) are the EGR valve
position (percent open) and VGT position (percent close). The control-oriented
airpath model also has an LPV form,

$\displaystyle z_{k+1}-z_{k+1}^{\tt
ss}(\rho_{k})=A(\rho_{k})\left[z_{k}-z_{k}^{\tt ss}(\rho_{k})\right]$
$\displaystyle+B(\rho_{k})\left[v_{k}-v_{k}^{\tt
ss}(\rho_{k})\right]+B_{f}(\rho_{k})\left[w_{{\tt inj},k}-w_{{\tt inj}}^{\tt
ss}(\rho_{k})\right],$ (6)

where $A,B:\mathbb{R}^{2}\to\mathbb{R}^{2\times 2}$, $z_{k}^{\tt
ss},v_{k}^{\tt ss}:\mathbb{R}^{2}\to\mathbb{R}^{2}$ are mappings that
determine equilibrium values of $z$ and $v$ corresponding to a given $\rho$,
$B_{f}:~{}\mathbb{R}^{2}\to\mathbb{R}^{2\times 1}$, and $w_{\tt inj}^{\tt
ss}(\rho)=\left[\begin{array}[]{cc}0&1\end{array}\right]\rho$. Note that the
model (2.3) includes $w_{{\tt inj},k}$ as an extra additive input. This input
has been added as includes it improves the model match in transients. Details
of the airpath LPV model identification and validation can be found in Zhang
et al. (2022a).

##  3 Integrated Emissions and Airpath MPC Design

The overall architecture of the proposed integrated emissions and airpath
controller is shown in Figure 6. For each engine speed ($N_{\tt e}$) and fuel
injection rate target ($w_{\tt inj}$), the outer-loop EMPC controller
generates adjusted target values for intake manifold pressure ($p_{\tt
im}^{adj}$), EGR rate ($\chi_{\tt egr}^{adj}$), as well as the fuel injection
command ($w_{\tt inj}^{adj}$) to reduce $NOx$ and enforce the constraint on
$Soot$. The fuel injection rate is applied directly to the engine while the
adjusted intake manifold pressure and EGR rate target values are passed to the
airpath MPC which tracks these targets by adjusting the EGR valve and VGT
positions, see Zhang et al. (2022a) for more details.

Figure 6: The block diagram of the integrated emissions and airpath control
system.

In the baseline engine control strategy, the intake manifold pressure and EGR
rate targets are static functions of the operating condition, i.e., $p_{\tt
im}^{\tt trg}(\rho):\mathbb{R}^{2}\to\mathbb{R}$ and $\chi_{\tt egr}^{\tt
trg}(\rho):\mathbb{R}^{2}\to\mathbb{R}$ that are implemented using look-up
tables. These targets are obtained during engine development and are chosen to
be optimal in steady-state operating conditions. In this paper, we focus on
the optimization of transient engine response, i.e., using MPC to shape the
transient response of the engine as it transitions between operating points.

###  3.1 EMPC for Emissions Control

Our economic MPC (EMPC) exploits a rate-based prediction model to enhance
closed-loop robustness and disturbance rejection. Rate-based MPC was applied
to aircraft turbofan engines by DeCastro (2007) and to diesel engine airpath
control by Huang et al. (2013, 2014, 2016).

First, in reference to the model (5), we define rate variables $\Delta x_{k}$
and $\Delta u_{k}$ as

$\Delta x_{k}=\tilde{x}_{k}-\tilde{x}_{k-1},\quad\Delta
u_{k}=\tilde{u}_{k}-\tilde{u}_{k-1},$

that correspond to the increments in the state and control variables,
respectively. Then, assuming $\rho_{k}$ remains constant over the prediction
horizon, the prediction model (5) implies that

$\Delta x_{k+1}=\tilde{A}(\rho_{k})\Delta x_{k}+\tilde{B}(\rho_{k})\Delta
u_{k}.$ (7)

The constant $\rho_{k}$ assumption is adopted widely in the literature on
engine and powertrain control with MPC (Norouzi et al. (2021)).

The rate-based reformulation as given in (7) is particularly effective when
applied to an LPV model, as the steady-state values of states and controls in
(4), i.e., $x_{k}^{ss}(\rho_{k})$ and $u_{k}^{ss}(\rho_{k})$, under the
assumption of $\rho_{k}$ remaining constant over the prediction horizon, does
not need to be known.

We let $\tilde{x}_{j|k}$ and $\tilde{u}_{j|k}$ denote the predicted state and
control values at time step $j,$ $0\leq j\leq N$, over the prediction horizon
when the prediction is made at the time step $k$. Then to be able to impose
constraints on $\tilde{x}_{j|k}$ and $\tilde{u}_{j|k}$ and compute the MPC
cost that penalizes the rate of change of input, $\Delta u_{k}$, we define the
augmented state vector,

$x_{j|k}^{\tt ext}=\left[\begin{array}[]{cccc}\Delta x_{j|k}^{\sf
T},&\tilde{x}_{j-1|k}^{\sf T},&\tilde{u}_{j-1|k}^{\sf
T}\end{array}\right]^{\sf T}.$

The rate-based prediction model (7) then implies

$\displaystyle x_{j+1|k}^{\tt ext}=\begin{bmatrix}\tilde{A}(\rho_{k})&0&0\\\
\mathbb{I}_{n_{x}\times n_{x}}&\mathbb{I}_{n_{x}\times n_{x}}&0\\\
0&0&\mathbb{I}_{n_{u}\times n_{u}}\end{bmatrix}x_{j|k}^{\tt ext}$
$\displaystyle\qquad\qquad+\begin{bmatrix}\tilde{B}(\rho_{k})\\\ 0\\\
\mathbb{I}_{n_{u}\times n_{u}}\end{bmatrix}\Delta u_{j|k},$ (8)

where $\Delta u_{j|k}$ is the control input in the extended system.

The EMPC is designed based on considerations of safety and drivability, so
that the fuel injection rate satisfies the specified upper and lower bounds,
as well as on the considerations of the emission requirements such as $Soot$
limits. The following optimal control problem is solved at each sampling
instant with the intake manifold pressure, EGR rate, and fuel injection rate
being optimized:

$\min_{\Delta u_{0|k},...,\Delta u_{N-1|k},\epsilon_{k}}J(\Delta
u,\epsilon,\rho)=\sum_{j=0}^{N}l(\Delta u_{j|k},\epsilon_{j|k},\rho_{k})$ (9a)

subject to

$\displaystyle\Delta x_{j+1|k}=\tilde{A}(\rho_{k})\Delta
x_{j|k}+\tilde{B}(\rho_{k})\Delta u_{j|k}$ (9b) $\displaystyle
u_{j|k}=u_{j-1|k}+\Delta u_{j|k},j={0},\ldots,N-1,$ (9c)
$\displaystyle{Soot}_{j|k}\leq{Soot}_{max}+\epsilon_{j|k},j=1,\ldots,N,$ (9d)
$\displaystyle 0.9w_{\tt inj}^{trg}\leq{w_{\tt inj}^{adj}}_{j|k}\leq w_{\tt
inj}^{trg},j=0,\ldots,N-1,$ (9e) $\displaystyle\epsilon_{j|k}\geq
0,j=0,\ldots,N-1,$ (9f)

where $N$ is the prediction horizon and $\epsilon_{k}$ is the slack variable
introduced to avoid infeasibility of the $Soot$ constraints. The stage cost
function $l$ is defined as

$\displaystyle l(\Delta u,\epsilon,\rho)=$ $\displaystyle\alpha(p_{\tt
im}^{trg}(\rho)-p_{\tt im}^{adj})^{2}+\beta(\chi_{\tt
egr}^{trg}(\rho)-\chi_{\tt egr}^{adj})^{2}+$ $\displaystyle\gamma(w_{\tt
inj}^{trg}-w_{\tt inj}^{adj})+\eta{NOx}+\zeta\epsilon+\Delta u^{\sf T}R\Delta
u$

where $\alpha,\beta,\gamma,\eta,\zeta>0$ and $R>0$ are tuning parameters and
reflect tracking objectives for the intake manifold pressure target, EGR rate
target, and fuel injection rate, a penalty on $NOx$ value, a penalty to soften
the $Soot$ constraint to guarantee feasibility, and a damping term. The fuel
tracking, $NOx$ term, and slack variable use 1-norm penalties since they are
more robust to ill-conditioning compared to quadratic penalties, facilitating
the solution of (9) numerically.

###  3.2 Tracking MPC for Airpath System

Our airpath MPC includes a feedforward (FF) MPC, to provide fast transient
response, and a feedback (FB) MPC that provides integral action and ensures
offset-free tracking. The FB MPC uses a rate-based formulation to eliminate
the steady-state error while the FF MPC only closes the loop around the model
of the plant. Details of design and tuning of the airpath MPC are presented in
Zhang et al. (2022a).

##  4 Simulation Results and Discussions

Four different tuning scenarios are considered to investigate how the
emissions can be altered according to different control requirements:

  * •

EMPC-A: Low $NOx$ penalty without $Soot$ limit

  * •

EMPC-B: High $NOx$ penalty without $Soot$ limit

  * •

EMPC-C: Low $NOx$ penalty with $Soot$ limit

  * •

EMPC-D: High $NOx$ penalty with $Soot$ limit

A baseline scenario without EMPC (w/o EMPC) was also generated, in which the
adjusted targets for $p_{\tt im}$, $\chi_{\tt egr}$ are the same as the values
from look-up table while $w_{\tt inj}$ is the same as defined by the operating
condition.

###  4.1 Case Study with Steps/Ramps

Figure 7 shows the simulation results at one operating point with different
$NOx$ penalties in the cost function during first, a fuel tip-in and tip-out
and then, a ramp change of engine speed. Here no $Soot$ limit was added. The
simulations are run for 30 $s$ first to achieve close to steady-state
temperature conditions before steps are applied. According to the results
summarized in Table 2, the EMPC reduces both cumulative $NOx$ and peak $NOx$
values. Figure 8 shows the corresponding intake manifold pressure, EGR rate,
and fuel injection rate values from EMPC-B simulation. The EMPC adjusts the
target $p_{\tt im}$ and $\chi_{\tt egr}$ to help reduce $NOx$ while the
airpath MPC is able to track these adjusted target values without steady-state
errors.

Figure 7: Comparison of (a) $NOx$ and (b) $Soot$ with different $NOx$
penalties during first, a fuel ($w_{\tt inj}$) tip-in and tip-out and then, a
ramp change of $N_{\tt e}$.  Figure 8: Comparison of actual , target from
look-up table, and adjusted target by EMPC (a) intake manifold pressure, (b)
EGR rate and (c) fuel injection rate with EMPC-B.

After adding the $Soot$ constraint, the EMPC is able to keep $Soot$ below a
predefined limit. Figure 9 shows the simulation results at the same operating
point and with the $Soot$ limit, and with different $NOx$ penalties in the
cost function. Table 3 indicates that both EMPC-C and EMPC-D have reduced
$Soot$ peak values. However, EMPC-D has higher average $Soot$ compared with
the baseline due to its large penalty on $NOx$ in the cost function. According
to Figure 10, $Soot$ values of both EMPCs are below the limit for most of the
time. Oscillations in the adjusted target values of $p_{\tt im}$ and
$\chi_{\tt egr}$ are observed in Figure 10 attributed to the model mismatch
between LPV and NN emissions models. The range of the oscillation in $p_{\tt
im}$ and $\chi_{\tt egr}$ is less than 0.01 bar and 1%, respectively, and is
relatively small. Such oscillation can also be reduced by using a larger
penalty $R$ on the damping term in the MPC cost function. At the same time,
the airpath MPC successfully tracks the adjusted target values.

Figure 9: Comparison of (a) $NOx$ and (b) $Soot$ with the same $NOx$ penalty and operating conditions as in Figure 7 with the additional of $Soot$ limit. Figure 10: Comparison of actual, target from look-up table, and adjusted target by EMPC (a) intake manifold pressure, (b) EGR rate, and (c) fuel injection rate with EMPC-D. Table 2: Comparison of the $NOx$ results for different EMPC scenarios. MPC |  Cumulative $NOx$ [%] |  Peak $NOx$ [%]
---|---|---
w/o EMPC | reference | reference
(baseline) |  |
EMPC-A |  $\downarrow$ 0.780% |  $\downarrow$ 5.700%
EMPC-B |  $\downarrow$ 6.695% |  $\downarrow$ 9.065%
EMPC-C |  $\uparrow$ 0.395% |  $\downarrow$ 5.341%
EMPC-D |  $\downarrow$ 3.033% |  $\downarrow$ 8.301%
Table 3: Comparison of the $Soot$ results for different EMPC scenarios. MPC |  Average $Soot$ [%] |  Peak $Soot$ [%]
---|---|---
w/o EMPC | reference | reference
(baseline) |  |
EMPC-A |  $\uparrow$ 1.067% |  $\uparrow$ 0.090%
EMPC-B |  $\uparrow$ 14.071% |  $\uparrow$ 14.946%
EMPC-C |  $\downarrow$ 1.300% |  $\downarrow$ 6.183%
EMPC-D |  $\uparrow$ 5.003% |  $\downarrow$ 4.765%

###  4.2 Transient Drive Cycle Simulations

We now evaluate the integrated emissions and airpath control strategy in more
complex simulation scenarios over the Federal Test Procedure (FTP) and WHTC.
The results are summarized in Tables 4 and 5. Both cycles run for 100 $s$
first to achieve close to steady-state temperature conditions before cycle
inputs are applied. Figures 11 and 12 show the comparison between EMPC-B and
EMPC-C against the baseline case over FTP and WHTC cycles, respectively. With
EMPC-B, $NOx$ are considerably reduced over both drive cycles with a sacrifice
of around 20% increase in average $Soot$. If the $Soot$ limit is imposed,
EMPC-C could reduce both average and peak $Soot$ in both drive cycles compared
with the baseline case. Moreover, according to Figures 11-(b) and 12-(b), most
of the peaks above the $Soot$ constraint are eliminated, illustrating EMPC-C
is able to reduce the visible smoke. Results of EMPC-D are only shown in
Tables 4 and 5, from which it can be seen that both cumulative $NOx$ and
average $Soot$ are reduced if the EMPC is made very aggressive on both $NOx$
reduction and $Soot$ limits.

Figure 11: Comparison of (a) $NOx$ and (b) $Soot$ control results over the FTP cycle. Figure 12: Comparison of (a) $NOx$ and (b) $Soot$ control results over the WHTC cycle. Table 4: Comparison of $NOx$ control results over FTP and WHTC cycles. MPC |  Cumulative $NOx$ [%] |  Peak $NOx$ [%]
---|---|---
w/o EMPC | reference | reference
(baseline) | FTP/WHTC | FTP/WHTC
EMPC-B |  $\downarrow$ 9.527%/$\downarrow$ 38.644% |  $\uparrow$ 8.374%/$\downarrow$ 33.530%
EMPC-C |  $\uparrow$ 4.056%/$\uparrow$ 2.963% |  $\uparrow$ 6.6059%/$\downarrow$ 1.649%
EMPC-D |  $\downarrow$ 0.451%/$\downarrow$ 15.378% |  $\uparrow$ 10.822%/$\downarrow$ 11.615%
Table 5: Comparison of $Soot$ control results over FTP and WHTC cycles. MPC |  Average $Soot$ [%] |  Peak $Soot$ [%]
---|---|---
w/o EMPC | reference | reference
(baseline) | FTP/WHTC | FTP/WHTC
EMPC-B |  $\uparrow$ 18.789%/$\uparrow$ 21.329% |  $\uparrow$ 17.800%/$\downarrow$ 7.238%
EMPC-C |  $\downarrow$ 14.268%/$\downarrow$ 9.521% |  $\downarrow$ 9.485%/$\downarrow$ 33.786%
EMPC-D |  $\downarrow$ 3.262%/$\downarrow$ 5.656% |  $\downarrow$ 9.147%/$\downarrow$ 34.741%

##  5 Conclusions

In this paper, a multi-layer Neural Network (NN) was trained and validated to
simultaneously predict $NOx$ and $Soot$ emissions of a diesel engine for both
transient and steady-state operating conditions. A co-simulation framework was
then implemented to integrate the NN model with a high-fidelity model of the
diesel engine airpath system built in GT-Power. An integrated control strategy
for emissions regulation and airpath system setpoint tracking was developed
based on model predictive control (MPC). The emission controller is based on
an economic MPC formulation that leverages a linear parameter varying (LPV)
model of $NOx$ and $Soot$ to minimize $NOx$ and enforce the $Soot$ limit. The
airpath controller is a rate-based MPC consisting of feedback and feedforward
loops, aiming to track adjusted intake manifold pressure and EGR rate
setpoints computed by the economic MPC.

The performance of our integrated emissions and airpath MPCs was assessed
through simulations over different steady-state and transient drive cycles in
terms of (i) reducing $NOx$, (ii) enforcing $Soot$ limit, and (iii) tracking
the adjusted intake manifold pressure and EGR rate setpoints. The results
highlight the potential of economic supervisory MPC to achieve improved
transient control and for reducing engine feedgas emissions. Future works
include using region-based EMPC calibration to improve transient performance
and replacing the LPV model used for emission prediction with a dynamic neural
network model to reduce model mismatch and improve the overall performance.

## References

  * Arrègle et al. (2008) Arrègle, J., López, J.J., Guardiola, C., and Monin, C. (2008).  Sensitivity study of a NOx estimation model for on-board applications.  Technical report, SAE Technical Paper.
  * Broomhead et al. (2016) Broomhead, T., Manzie, C., Hield, P., Shekhar, R., and Brear, M. (2016).  Economic model predictive control and applications for diesel generators.  _IEEE Transactions on Control Systems Technology_ , 25(2), 388–400.
  * De Cesare and Covassin (2011) De Cesare, M. and Covassin, F. (2011).  Neural network based models for virtual no x sensing of compression ignition engines.  Technical report, SAE Technical Paper.
  * DeCastro (2007) DeCastro, J.A. (2007).  Rate-based model predictive control of turbofan engine clearance.  _Journal of Propulsion and Power_ , 23(4), 804–813.
  * Hagena et al. (2006) Hagena, J.R., Filipi, Z.S., and Assanis, D.N. (2006).  _Transient diesel emissions: analysis of engine operation during a tip-in_.  SAE International New York, NY, USA.
  * Heywood (2018) Heywood, J.B. (2018).  _Internal combustion engine fundamentals_.  McGraw-Hill Education.
  * Huang et al. (2020) Huang, C., Salehi, R., Ersal, T., and Stefanopoulou, A.G. (2020).  An energy and emission conscious adaptive cruise controller for a connected automated diesel truck.  _Vehicle System Dynamics_ , 58(5), 805–825.
  * Huang et al. (2018) Huang, M., Liao-McPherson, D., Kim, S., Butts, K., and Kolmanovsky, I. (2018).  Toward real-time automotive model predictive control: A perspective from a diesel air path control development.  In _2018 Annual American Control Conference_ , 2425–2430. IEEE.
  * Huang et al. (2014) Huang, M., Nakada, H., Butts, K., and Kolmanovsky, I. (2014).  Robust rate-based model predictive control of diesel engine air path.  In _American Control Conference_ , 1505–1510.
  * Huang et al. (2013) Huang, M., Nakada, H., Polavarapu, S., Butts, K., and Kolmanovsky, I. (2013).  Rate-based model predictive control of diesel engines.  _IFAC Proceedings Volumes_ , 46(21), 177–182.
  * Huang et al. (2016) Huang, M., Zaseck, K., Butts, K., and Kolmanovsky, I. (2016).  Rate-based model predictive controller for diesel engine air path: Design and experimental evaluation.  _IEEE Transactions on Control Systems Technology_ , 24(6), 1922–1935.
  * Li et al. (2017) Li, H., Butts, K., Zaseck, K., Liao-McPherson, D., and Kolmanovsky, I. (2017).  Emissions modeling of a light-duty diesel engine for model-based control design using multi-layer perceptron neural networks.  Technical report, SAE Technical Paper.
  * Liao-McPherson et al. (2020) Liao-McPherson, D., Huang, M., Kim, S., Shimada, M., Butts, K., and Kolmanovsky, I. (2020).  Model predictive emissions control of a diesel engine airpath: Design and experimental evaluation.  _International Journal of Robust and Nonlinear Control_ , 30(17), 7446–7477.
  * Liu et al. (2021) Liu, Z., Dizqah, A.M., Herreros, J.M., Schaub, J., and Haas, O. (2021).  Simultaneous control of nox, soot and fuel economy of a diesel engine with dual-loop egr and vnt using economic mpc.  _Control Engineering Practice_ , 108, 104701.
  * Norouzi et al. (2021) Norouzi, A., Heidarifar, H., Shahbakhti, M., Koch, C.R., and Borhan, H. (2021).  Model predictive control of internal combustion engines: A review and future directions.  _Energies_ , 14(19), 6251.
  * Ortner and Del Re (2007) Ortner, P. and Del Re, L. (2007).  Predictive control of a diesel engine air path.  _IEEE Transactions on Control Systems Technology_ , 15(3), 449–456.
  * Paszke et al. (2019) Paszke, A., Gross, S., Massa, F., Lerer, A., Bradbury, J., Chanan, G., Killeen, T., Lin, Z., Gimelshein, N., Antiga, L., et al. (2019).  Pytorch: An imperative style, high-performance deep learning library.  _Advances in Neural Information Processing Systems_ , 32.
  * Prokhorov (2008) Prokhorov, D. (2008).  Neural networks in automotive applications.  In _Computational intelligence in automotive applications_ , 101–123. Springer.
  * Winkler-Ebner et al. (2010) Winkler-Ebner, B., Hirsch, M., Del Re, L., Klinger, H., and Mistelberger, W. (2010).  Comparison of virtual and physical nox-sensors for heavy duty diesel engine application.  _SAE International Journal of Engines_ , 3(1), 1124–1139.
  * Zhang et al. (2022a) Zhang, J., Amini, M.R., Kolmanovsky, I., Tsutsumi, M., and Nakada, H. (2022a).  Benefits of feedforward for model predictive airpath control of diesel engines.  In _10th IFAC Symposium on Robust Control Design_ , 181–186.
  * Zhang et al. (2022b) Zhang, J., Amini, M.R., Kolmanovsky, I., Tsutsumi, M., and Nakada, H. (2022b).  Development of a model predictive airpath controller for a diesel engine on a high-fidelity engine model with transient thermal dynamics.  In _American Control Conference_ , 3058–3063.


#  Collectives: Towards Cloud-aware Collectives for ML Workloads with Rank
Reordering

Liang Luo∗‡, Jacob Nelson†, Arvind Krishnamurthy∗, Luis Ceze∗
∗University of Washington, †Microsoft Research, ‡Facebook AI

###### Abstract

ML workloads are becoming increasingly popular in the cloud. Good cloud
training performance is contingent on efficient parameter exchange among VMs.
We find that Collectives, the widely used distributed communication
algorithms, cannot perform optimally out of the box due to the hierarchical
topology of datacenter networks and multi-tenancy nature of the cloud
environment.

In this paper, we present Collectives (Cloud-aware collectives), a prototype
that accelerates collectives by reordering the ranks of participating VMs such
that the communication pattern dictated by the selected collectives operation
best exploits the locality in the network. Collectives is non-intrusive,
requires no code changes nor rebuild of an existing application, and runs
without support from cloud providers. Our preliminary application of
Collectives on allreduce operations in public clouds results in a speedup of
up to 3.7x in multiple microbenchmarks and 1.3x in real-world workloads of
distributed training of deep neural networks and gradient boosted decision
trees using state-of-the-art frameworks.

##  I Introduction

Collective communication operations are an essential component of many data-
parallel computation frameworks. Originally developed for high-performance
computing frameworks such as MPI [49], they are now widely used for cloud-
based distributed machine learning [48, 9, 1, 33, 27, 8, 43] workloads. With
increasingly more complex models [37, 35, 42] calling for larger data
exchanges, and rapid deployment of faster accelerators [11, 21, 26] demanding
more frequent exchanges, efficient execution of these workloads contingent on
efficient collectives.

Unfortunately, achieving good collectives performance in a cloud environment
is fundamentally more challenging than in an HPC world, because the user has
no control over node placement, topology and has to share the infrastructure
with other tenants. These constraints have a strong implication on the
performance of collectives. As a result, the bottleneck of running these
workloads on the cloud has shifted from computation to communication [50, 30,
29].

Figure 1: Performance distribution of allreduce task of 100MB data with ring
algorithm varies widely with 500 random rank orders on Azure. Figure 2:
Pairwise RTT probe on 64 Azure Standard F64 v2 VMs shows non-uniform latency
in point to point VM communication.

Consider a common practice of applying allreduce ring collectives, a popular
algorithm, in the cloud context, where a randomly-ordered IP list (obtained
through the provider) of VMs is used to form a virtual ring on which data is
passed along, with $i$-th VM sending data to $i+1$-th VM. But do different
ways of forming ring (through permutation of VMs in the list) exhibit the same
performance? The answer is most likely no, as the ring that corresponds to
shorter total hop cost will likely perform better (Figure 1). On the other
hand, not all ways of forming rings achieve the same cost, because the point
to point communication cost (bandwidth, latency, or collectively referred to
as locality in this work) is different across VMs (Figure 2), due to the
hierarchical structure of the datacenter network, and the dynamic nature of
traffic from other tenants. Consequently, running collectives with a randomly-
ordered list of VMs results in unpredictable and subpar performance.

Our work focuses on discovering a permutation of the IP list that exploits the
network locality for efficient communication, in a completely transparent way,
by minimizing the cost model of a given collectives parameterized with the
actual hop cost. To do so, we need to (1) efficiently identify the underlying
network constraints (or collectively, locality); (2) accurately build cost
model for the collectives at hand; (3) effectively approximate the minimum of
complex cost functions.

This paper proposes Collectives, a tool that uses network probes to discover
locality within the underlying datacenter network, and uses it to solve a
communication cost minimization problem with constraints, with the rank of
each VM as the unknowns. We use reordered ranks as input to unmodified
communication backends in microbenchmarks including OMB [17], Nvidia NCCL [6],
Facebook Gloo [23] and real-world workloads of training deep neural networks
with Pytorch/Caffe2 and gradient boosted decision trees using LightGBM [33,
27] and our preliminary results show a speedup of up to 3.7x in various
allreduce operations and 1.3x in end-to-end performance across EC2 and Azure.

##  II Background

We provide an overview of typical structures and performance implications for
datacenter networks and a brief introduction to the various popular
collectives operations.

###  II-A Locality in Datacenter Network

Modern datacenters are hierarchical, with machines connecting to a top of rack
switch, which are in turn connected to upper-level devices [36, 24, 44, 28].
This particular topology induces locality [30], as the communication
performance between two physical hosts is not the same. For example, VMs
within the same rack have the best and stable performance, as the physical
link is not shared. On the other hand, links between hosts residing in
different racks are shared, and the communication performance depends on
factors like hop count, link congestion, oversubscription ratio [15], and
dynamic load. Topology information is crucial for achieving optimal
performance as many collectives implementations generate routines based on
this information [19, 38, 39, 20, 29]. But in a cloud-environment, this
information is hidden. Various attempts are made to reconstruct the physical
affinity, e.g., PLink [30] uses DPDK-based latency probes and K-Means
clustering to find hosts with high physical affinity.

###  II-B Collectives

Collectives works by decomposing an operation into a series of point to point
operations between two nodes according to a predefined communication schedule.
collectives most often appear in MPI contexts [45, 47, 41, 16, 13]. Typically,
all nodes in collectives participate in the communication, usually running
symmetric tasks. collectives can be used in many tasks such as (all)reduce,
broadcast, (all)gather, and (reduce)scatter [7], and it is thus impractical to
individually optimize each task. Fortunately, many of the tasks can be
decomposed into multiple stages of collectives primitives (e.g., allreduce can
be decomposed into reducescatter followed by allgather. Therefore, we only
need to focus on accelerating such primitives. We now introduce these
algorithms. We use $N$ as the number of participating nodes, $S$ as the amount
of data to process per node.

Figure 3: Rounds of communications (color-coded) in various popular algorithms
used in collectives with 4 nodes.

Ring [40]. As shown in Figure 3(a), ring algorithms work by connecting nodes
to form a virtual ring. Data is then passed along the ring sequentially. Ring
algorithms require $O(N)$ steps to complete, sending $O(NS)$ amount of data.

Halving Doubling [47]. As shown in Figure 3(b), halving doubling works by
recursively doubling the distance (in terms of rank ID) while halving the
total amount of data sent in each round, requiring $O(log_{2}{N})$ steps to
finish while sending $O(NS)$ amount of data on the wire.

Tree. In one form, a single tree is built where data is transferred from
leaves to the root and vice-versa [25]; in a more optimized setting, a pair of
complementary binary trees are built to fully utilize the full bisection
bandwidth [46], each sending and receiving $S/2$. For binary trees,
$O(log_{2}N)$ rounds of communication are required, sending $O(NS)$ bytes.

BCube [3]. BCube is very similar to halving doubling from a structural
perspective, in the sense that nodes are organized into a group of $B$ peers.
BCube operates in $O(log_{B}N)$ rounds, and each node in each round would peer
with a unique node in another $B-1$ groups. Each node communicates
$\frac{S}{B^{i}}$ amount of data in round $i$. BCube achieving a total bytes
on wire of $O(\sum_{i=0}^{log_{B}N-1}\frac{S}{B^{i+1}})$.

##  III Motivation

This section motivates Collectives by demonstrating the implication of rank
order on the performance of cloud-based collectives. We start by highlighting
the asymmetric, non-uniform link cost in the cloud environment, by launching
64 Standard F64 VMs on the Azure cloud. We then run an in-house hybrid DPDK
and ping-based [4] latency probe (§IV-B) between each pair of VM node, using
the same technique in PLink. The result is summarized in Figure 2 as a
heatmap. We observe the pairwise round-trip latency can range from sub-10 to
hundreds of microseconds.

We proceed to examine the performance of allreduce operation using Gloo [2],
running the Ring chunked algorithm with 512 Standard F16 VMs. To derive a
performance distribution, we use 500 randomly generated rank orders to
generate 500 samples, and each is the average runtime of a 10-iteration
reduction of 100MB of data. The result is summarized in the yellow
distribution in Figure 1. The performance of different rank orderings of VM
varies drastically, ranging from 330ms to 3400ms, with a mean of 1012ms and a
standard deviation of 418ms. Now we have established the profound influence of
rank ordering of VM nodes on the performance of collectives algorithms, the
goal of this paper is to derive an approximately optimal rank-ordering given a
selected collectives algorithm such that it maximizes the performance.

##  IV Design and Implementation

We now describe Collectives, a tool that takes in a list of VM nodes and a
target algorithm, accurately and efficiently probes their pairwise distance,
and uses that information to construct a rank order of VMs that attempts to
minimize the total cost of communication.

###  IV-A Cost Models for Collective Algorithms

Collectives builds a cost model $\mathbb{C}_{\mathbb{O}}$ for each popular
algorithm used in collectives $\mathbb{O}$, parameterized with the number of
participating nodes $N$ and size $S$. This section details the cost models for
popular algorithms. We use $c_{i,j}(S)$ to refer to the cost for transferring
$S$ amount of data from node $i$ to $j$. We further define
$MAX_{i=0}^{j}(f(i))=MAX(f(0),...,f(j))$. We assume $N$ a power of 2 to
simplify explanation, and allow arbitrary rank $r$ to alias to canonical rank
$(r+N)\text{ mod }N$.

Ring. The cost model of the ring algorithm is the sum of the cost of each hop
when traversing the ring:

$\mathbb{C}_{r}(N,c,S)=\sum_{i=0}^{N-1}c_{i,i-1}(S)$

Having Doubling. The cost of halving doubling is the sum of costs for each
round of communication, which in turn is the max cost of all communications in
that round.

$\mathbb{C}_{hd}(N,c,S)=\sum_{i=0}^{log_{2}N-1}MAX_{j=0}^{\frac{N}{2}-1}c_{j,j+2^{i}}(\frac{S}{2^{i+1}})$

Tree. The cost of running tree algorithms depends on the number of trees and
how trees are constructed. The total cost is the maximum cost of all trees,
which is in turn determined by the maximum cost of each subtree. We provide a
cost model for a popular variant of tree algorithm: double binary tree as used
in [5].

$\mathbb{C}_{dbt}(N,c,S)=T(0,N-1,S)$

where $T(i,j,S)$ is expressed recursively:

$\displaystyle T(i,j,S)=\left\\{\begin{array}[]{l@{}l}0\textbf{ if $i\geq
j$}\\\
MAX(c_{\frac{i+j}{2},\frac{3i+j}{2}-1}(\frac{S}{2})+T(i,\frac{i+j}{2}-1),\\\
c_{\frac{i+j}{2},\frac{i+3j}{2}+1}(\frac{S}{2})+T(\frac{i+j}{2}+1,j))\textbf{
otherwise}\end{array}\right.$

Similarly a mirrored tree is built by decrementing each node’s rank in the
tree without changing the tree structure.

BCube. The cost of running the BCube algorithm is similar to halving doubling,
except in each round, each node communicates with $B-1$ peers, instead of 1.

$\mathbb{C}_{b}(N,c,S,B)=\sum_{i=0}^{log_{B}N-1}MAX_{j=0}^{\frac{N}{B}-1}MAX_{k=1}^{B}c_{j,j+kB^{i}}(\frac{S}{B^{i+1}})$

###  IV-B Probing for Pairwise Distance

We need to determine values for $c_{i,j}(S)$ with end-to-end measurements. In
this work, we use a latency-centric view for the cost component. The rationale
behind this stems from the well-known theoretical TCP bandwidth model of
$BW=O(\frac{MSS}{RTT\sqrt{p}})$ [32]: given constant drop rate $p$ and window
$MSS$, higher latency induces lower bandwidth in TCP streams. This
conveniently lets us approximate costs by only probing for latency. We adopt
the probing pipeline used in PLink, which focuses on discovering physical
locality with an in-house DPDK based echo tool, leveraging network enhancement
provided by the clouds [34, 12]. Each pair of nodes receive a total of $10k$
probes from sequentially and bidirectionally. To derive an accurate reading,
we take the RTT of 10th percentile to filter out interference during probes.
Each probe is a UDP packet with a 32-bit payload that encodes sequence number
and round id for fault tolerance. When DPDK cannot be used, we use fping, a
ICMP Echo-based latency probing tool. For each entry in $c$, we update
$c_{i,j}\leftarrow MAX(c_{i,j},c_{j,i})$ to make it symmetric.

###  IV-C Minimizing the Cost Model

We parameterize the cost model with values of probed $c$. To derive a rank
ordering that minimizes $\mathbb{C_{O}}$, we perform the following
transformation: let set of variables $\mathbb{R}$ defined as
$r_{i},i\in[0,N-1]$ be a permutation of $[0,N-1]$ to be solved, and we replace
each $c_{i,j}$ with $c_{r_{i},r_{j}}$. We can then establish a bijection from
the original rank ordering to the desired order $r_{i}\leftrightarrow i$ once
$r_{i}$s are solved. We flatten $c_{i,j}\leftrightarrow c^{\prime}_{iN+j}$ to
use theory of arrays to allow direct solving with conventional optimizing SMT
solvers such as Z3 [18, 22].

Unfortunately, we find solvers inefficient, perhaps due to the non-convex,
non-linear nature of the objective function and a large search space ($N!$).
Thus, we take a two-stage solving process. The first step employs a range of
stochastic search techniques such as simulated annealing [14], with a few
standard heuristics (e.g., permuting a random sub-array, permuting random
pairs) for obtaining neighboring states and a timeout. When the search returns
with an initial result $C_{0}$, we generate an additional SMT constraint
$\mathbb{C_{O}}<C_{0}$ to better guide pruning for solvers. We let the solver
continue to run for a few minutes, and we either find a better solution or
will use $C_{0}$ as the final value. The end-product of this process is a
rearranged list of VMs.

##  V Preliminary Evaluation

We evaluate Collectives with a series of microbenchmarks from various
communication backends and real-world applications that use collectives. We
represent speedup by comparing the performance we get from the best rank order
and the worst rank order. We avoid comparing with the original rank order
because it is random and unstable (Figure 1).

###  V-A Experimental Setup

Our experiments are conducted on two public clouds, Azure and EC2. We enable
network acceleration on both clouds and set TCP congestion control protocol to
DCTCP. We include microbenchmarks that exercise ring, having doubling, double
binary tree, and Bcube algorithms. All experiments run on Ubuntu 19. We focus
our evaluation on one of the most important tasks in collectives, allreduce
for its popularity.

###  V-B Prediction Accuracy of Cost Model

While the goal of the cost model is not to predict the actual performance, but
rather, it should preserve the relative order of performance, i.e.,
$p_{pred}(\mathbb{R}_{1})<p_{pred}(\mathbb{R}_{2})\implies
p_{real}(\mathbb{R}_{1})<p_{real}(\mathbb{R}_{2})$ should hold true for as
many pairs of $(R_{i},R_{j})$s as possible. We demonstrate this for ring based
collectives by generating 10 different rank orders, with the $i$-th order
approximately corresponds to the $10i$-th percentile in the range of costs
found by the solver. We obtain performance data for Facebook Gloo and OpenMPI
running OSU Benchmark on 64 F16 nodes on Azure and 64 C5 nodes on EC2. We then
compute Spearman [10] correlation coefficient between the predicted
performance and the actual performance for each setup (Table I). We found the
cost function predicting and actual collectives performance exhibits strong
correlation.

Setup | Azure | EC2
---|---|---
Gloo Ring 100MB | 0.58 | 0.78
OpenMPI Ring 100MB | 0.81 | 0.94
Table I: Spearman correlation coefficient between predicted performance from
cost model and actual performance.

###  V-C Microbenchmark Performance

We evaluate Collectives’s efficacy with microbenchmarks of algorithms
introduced in II. We report a mean speedup of 20 iterations. We run all
benchmarks with 512 F16 nodes on Azure, except for NCCL, which runs on 64
P3.8xLarge GPU nodes on EC2. Specifically, we set $B=4$ for BCube; for NCCL,
we use a single binary tree reduction for small buffers and a ring for large
buffers. In all benchmarks, we reduce a buffer of 100MB, except for Nvidia
NCCL, where we reduce a small buffer of 4B to trigger the tree algorithm.

Figure 4: Speedups achieved by just using the rank ordering produced by our
tool in various algorithms across multiple distributed communication backends
at a large scale.

Figure 4 shows a summary of speedups achieved using Collectives on these
benchmarks with the ring-family algorithms benefiting the most, up to 3.7x. We
speculate the reason for the effectiveness is that they have a much wider
performance distribution, as each permutation of the order can potentially
generate a different performance (cost of each hop is on the critical path);
they also have simpler cost model, allowing the solvers to quickly navigate
the objective landscape. On the other hand, halving doubling, BCube, and tree
algorithms have complex objectives – sum of maximums, resulting in a narrower
performance distribution because mutation of the ordering may not change the
cost at all if the mutation does not cause the critical path to change.

###  V-D End-to-end Performance Impact on Real-world Applications

Distributed Gradient Boosted Decision Tree. We evaluate Collectives’s impact
on LightGBM [27], a gradient boosted decision tree training system. We use
data parallelism to run lambdarank with metric ndcg. Communication-wise, this
workload runs two tasks: allreduce and reducescatter, and they are called
sequentially in each split of the iteration. At our scale of 512 nodes,
LightGBM automatically chooses to use halving and doubling for both
reducescatter and allreduce. We use a dataset that represents an actual
workload in our commercial setting with 5K columns and a total size of 10GB
for each node. We train 1K trees, each with 120 leaves. We exclude the time it
takes to load data from disk to memory and report an average speedup of 1000
iterations. Collectives generated rank orders speed up training by 1.3x.

Distributed Deep Neural Network. We show Collectives’s effectiveness on
distributed training of DNNs with Caffe2/Pytorch, on 64 EC2 p3.8xLarge nodes
with data parallelism and a batch size of 64/GPU. We train AlexNet on the
ImageNet dataset. Since our Collectives does not change computation and only
improves communication efficiency of the allreduce operation at iteration
boundary, we report speedup of training in terms of images/second, averaged
across 50 iterations. We use the ring chunked algorithm, which achieves the
best baseline performance, and the Collectives-optimized rank ordering of VM
nodes achieves a speedup of 1.2x.

##  VI Discussion and Limitations

Generalizability Study. Due to lack of resources, we only evaluated
Collectives on a few VM allocations on each cloud. Further investigation is
needed to understand how Collectives performs across different physical
allocations and whether Collectives’s rank order reacts to cloud variance
well.

Limitations to Cost Modeling. Our use of latency-centric cost function is
reasonable and effective, but is not perfect: it is possible that different
bandwidth can correspond to the same latency depending on link condition. This
may cause unoptimal solution when a transfer is bandwidth-bound. Further study
is needed to determine how to incorporate bandwidth related cost into the cost
function without incurring high probing overhead.

Complements to Cost Models. While building accurate cost models is difficult,
we can dynamically adjust rank ordering with help from tools such as TCP_INFO
[31] that monitors link properties such as latency and bandwidth. Since we
know the communication pattern, we can determine the critical path and find
bottleneck transfer between node $n_{i}$ and $n_{j}$ in the system. From
there, we can find a $n_{k}$ to replace $n_{i}$ such that the replacement
results in a minimized cost objective.

Adapting to Dynamic Traffic. The above mechanism can be applied to adapt to
dynamic network load in the cloud environment. The framework, however, must
support the dynamic change of node ranks, which is possible and should come
with a small cost as a full mesh of connections can be established beforehand
among all nodes.

## References

  * [1] “Extend mxnet distributed training by mpi allreduce - mxnet - apache software foundation,” https://cwiki.apache.org/confluence/display/MXNET/Extend+MXNet+Distributed+Training+by+MPI+AllReduce, (Accessed on 01/15/2020).
  * [2] “Github - facebookincubator/gloo: Collective communications library with various primitives for multi-machine training.” https://github.com/facebookincubator/gloo, (Accessed on 01/24/2020).
  * [3] “gloo/algorithms.md at master · facebookincubator/gloo,” https://github.com/facebookincubator/gloo/blob/master/docs/algorithms.md, (Accessed on 01/26/2020).
  * [4] “Home - dpdk,” https://www.dpdk.org/, (Accessed on 01/24/2020).
  * [5] “Massively scale your deep learning training with nccl 2.4 — nvidia developer blog,” https://devblogs.nvidia.com/massively-scale-deep-learning-training-nccl-2-4/, (Accessed on 01/25/2020).
  * [6] “Nvidia collective communications library (nccl) — nvidia developer,” https://developer.nvidia.com/nccl, (Accessed on 01/19/2020).
  * [7] “pages.tacc.utexas.edu/~eijkhout/pcse/html/mpi-collective.html,” http://pages.tacc.utexas.edu/~eijkhout/pcse/html/mpi-collective.html, (Accessed on 01/22/2020).
  * [8] “Spark-mpi: Approaching the fifth paradigm - databricks,” https://databricks.com/session/spark-mpi-approaching-the-fifth-paradigm, (Accessed on 01/15/2020).
  * [9] “Writing distributed applications with pytorch — pytorch tutorials 1.3.1 documentation,” https://pytorch.org/tutorials/intermediate/dist_tuto.html, (Accessed on 01/15/2020).
  * [10] _Spearman Rank Correlation Coefficient_. New York, NY: Springer New York, 2008, pp. 502–505. [Online]. Available: https://doi.org/10.1007/978-0-387-32833-1_379
  * [11] Amazon, “Amazon ec2 f1 instances,” https://aws.amazon.com/ec2/instance-types/f1/, (Accessed on 03/04/2019).
  * [12] Amazon, “Enable and configure enhanced networking for ec2 instances,” https://aws.amazon.com/premiumsupport/knowledge-center/enable-configure-enhanced-networking/, (Accessed on 01/06/2019).
  * [13] V. Bala, J. Bruck, R. Cypher, P. Elustondo, A. Ho, C.-T. Ho, S. Kipnis, and M. Snir, “Ccl: A portable and tunable collective communication library for scalable parallel computers,” _IEEE Transactions on Parallel and Distributed Systems_ , vol. 6, no. 2, pp. 154–164, 1995.
  * [14] D. Bertsimas and J. Tsitsiklis, “Simulated annealing,” _Statistical Science_ , vol. 8, 02 1993.
  * [15] K. Bilal, S. U. Khan, J. Kolodziej, L. Zhang, K. Hayat, S. A. Madani, N. Min-Allah, L. Wang, and D. Chen, “A comparative study of data center network architectures,” in _ECMS_ , 2012.
  * [16] E. K. Blum, X. Wang, and P. Leung, “Architectures and message-passing algorithms for cluster computing: Design and performance,” _Parallel Computing_ , vol. 26, no. 2-3, pp. 313–332, 2000.
  * [17] D. Bureddy, H. Wang, A. Venkatesh, S. Potluri, and D. K. Panda, “Omb-gpu: A micro-benchmark suite for evaluating mpi libraries on gpu clusters,” in _Recent Advances in the Message Passing Interface_ , J. L. Träff, S. Benkner, and J. J. Dongarra, Eds. Berlin, Heidelberg: Springer Berlin Heidelberg, 2012, pp. 110–120.
  * [18] L. De Moura and N. Bjørner, “Z3: An efficient smt solver,” in _International conference on Tools and Algorithms for the Construction and Analysis of Systems_. Springer, 2008, pp. 337–340.
  * [19] A. Faraj, P. Patarasuk, and X. Yuan, “Bandwidth efficient all-to-all broadcast on switched clusters,” in _2005 IEEE International Conference on Cluster Computing_ , Sep. 2005, pp. 1–10.
  * [20] A. Faraj and X. Yuan, “Message scheduling for all-to-all personalized communication on ethernet switched clusters,” in _19th IEEE International Parallel and Distributed Processing Symposium_ , April 2005, pp. 10 pp.–.
  * [21] J. Fowers, K. Ovtcharov, M. Papamichael, T. Massengill, M. Liu, D. Lo, S. Alkalay, M. Haselman, L. Adams, M. Ghandi _et al._ , “A configurable cloud-scale dnn processor for real-time ai,” in _Proceedings of the 45th Annual International Symposium on Computer Architecture_. IEEE Press, 2018, pp. 1–14.
  * [22] Google, “Or-tools — google developers,” https://developers.google.com/optimization, (Accessed on 01/26/2020).
  * [23] P. Goyal, P. Dollár, R. Girshick, P. Noordhuis, L. Wesolowski, A. Kyrola, A. Tulloch, Y. Jia, and K. He, “Accurate, large minibatch sgd: Training imagenet in 1 hour,” 2017.
  * [24] A. Greenberg, J. R. Hamilton, N. Jain, S. Kandula, , P. Lahiri, D. Maltz, and and, “Vl2: A scalable and flexible data center network.” Association for Computing Machinery, Inc., August 2009, recognized as one of ”the most important research results published in CS in recent years”. [Online]. Available: https://www.microsoft.com/en-us/research/publication/vl2-a-scalable-and-flexible-data-center-network/
  * [25] F. N. Iandola, M. W. Moskewicz, K. Ashraf, and K. Keutzer, “Firecaffe: near-linear acceleration of deep neural network training on compute clusters,” in _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition_ , 2016, pp. 2592–2600.
  * [26] N. P. Jouppi, C. Young, N. Patil, D. Patterson, G. Agrawal, R. Bajwa, S. Bates, S. Bhatia, N. Boden, A. Borchers, R. Boyle, P.-l. Cantin, C. Chao, C. Clark, J. Coriell, M. Daley, M. Dau, J. Dean, B. Gelb, T. V. Ghaemmaghami, R. Gottipati, W. Gulland, R. Hagmann, C. R. Ho, D. Hogberg, J. Hu, R. Hundt, D. Hurt, J. Ibarz, A. Jaffey, A. Jaworski, A. Kaplan, H. Khaitan, D. Killebrew, A. Koch, N. Kumar, S. Lacy, J. Laudon, J. Law, D. Le, C. Leary, Z. Liu, K. Lucke, A. Lundin, G. MacKean, A. Maggiore, M. Mahony, K. Miller, R. Nagarajan, R. Narayanaswami, R. Ni, K. Nix, T. Norrie, M. Omernick, N. Penukonda, A. Phelps, J. Ross, M. Ross, A. Salek, E. Samadiani, C. Severn, G. Sizikov, M. Snelham, J. Souter, D. Steinberg, A. Swing, M. Tan, G. Thorson, B. Tian, H. Toma, E. Tuttle, V. Vasudevan, R. Walter, W. Wang, E. Wilcox, and D. H. Yoon, “In-datacenter performance analysis of a tensor processing unit,” in _Proceedings of the 44th Annual International Symposium on Computer Architecture_ , ser. ISCA ’17. New York, NY, USA: ACM, 2017, pp. 1–12. [Online]. Available: http://doi.acm.org/10.1145/3079856.3080246
  * [27] G. Ke, Q. Meng, T. Finley, T. Wang, W. Chen, W. Ma, Q. Ye, and T.-Y. Liu, “Lightgbm: A highly efficient gradient boosting decision tree,” in _NIPS_ , 2017.
  * [28] M. Liu, L. Luo, J. Nelson, L. Ceze, A. Krishnamurthy, and K. Atreya, “IncBricks: Toward in-network computation with an in-network cache,” _SIGOPS Oper. Syst. Rev._ , vol. 51, no. 2, pp. 795–809, Apr. 2017. [Online]. Available: http://doi.acm.org/10.1145/3093315.3037731
  * [29] L. Luo, J. Nelson, L. Ceze, A. Phanishayee, and A. Krishnamurthy, “Parameter hub: A rack-scale parameter server for distributed deep neural network training,” in _Proceedings of the ACM Symposium on Cloud Computing_ , ser. SoCC ’18. New York, NY, USA: Association for Computing Machinery, 2018, p. 41–54. [Online]. Available: https://doi.org/10.1145/3267809.3267840
  * [30] L. Luo, P. West, J. Nelson, A. Krishnamurthy, and L. Ceze, “Plink: Discovering and exploiting locality for accelerated distributed training on the public cloud,” in _MLSys 2020_ , 2020.
  * [31] M. Mathis, J. Heffner, and R. Reddy, “Web100: extended tcp instrumentation for research, education and diagnosis,” _ACM SIGCOMM Computer Communication Review_ , vol. 33, no. 3, pp. 69–79, 2003.
  * [32] M. Mathis, J. Semke, J. Mahdavi, and T. Ott, “The macroscopic behavior of the tcp congestion avoidance algorithm,” _ACM SIGCOMM Computer Communication Review_ , vol. 27, no. 3, pp. 67–82, 1997.
  * [33] Q. Meng, G. Ke, T. Wang, W. Chen, Q. Ye, Z.-M. Ma, and T.-Y. Liu, “A communication-efficient parallel algorithm for decision tree,” in _Advances in Neural Information Processing Systems 29_ , D. D. Lee, M. Sugiyama, U. V. Luxburg, I. Guyon, and R. Garnett, Eds. Curran Associates, Inc., 2016, pp. 1279–1287. [Online]. Available: http://papers.nips.cc/paper/6381-a-communication-efficient-parallel-algorithm-for-decision-tree.pdf
  * [34] Microsoft, “Create an azure virtual machine with accelerated networking — microsoft docs,” https://docs.microsoft.com/en-us/azure/virtual-network/create-vm-accelerated-networking-cli, (Accessed on 01/06/2019).
  * [35] Microsoft, “Turing-nlg: A 17-billion-parameter language model by microsoft - microsoft research,” https://www.microsoft.com/en-us/research/blog/turing-nlg-a-17-billion-parameter-language-model-by-microsoft/, (Accessed on 05/25/2020).
  * [36] R. N. Mysore, A. Pamboris, N. Farrington, N. Huang, P. Miri, S. Radhakrishnan, V. Subramanya, and A. Vahdat, “Portland: a scalable fault-tolerant layer 2 data center network fabric,” in _SIGCOMM_ , 2009.
  * [37] OpenAI, “Ai and compute,” https://openai.com/blog/ai-and-compute/, (Accessed on 05/25/2020).
  * [38] P. Patarasuk, A. Faraj, and Xin Yuan, “Pipelined broadcast on ethernet switched clusters,” in _Proceedings 20th IEEE International Parallel Distributed Processing Symposium_ , April 2006, pp. 10 pp.–.
  * [39] P. Patarasuk and X. Yuan, “Bandwidth efficient all-reduce operation on tree topologies,” in _2007 IEEE International Parallel and Distributed Processing Symposium_ , March 2007, pp. 1–8.
  * [40] P. Patarasuk and X. Yuan, “Bandwidth optimal all-reduce algorithms for clusters of workstations,” _Journal of Parallel and Distributed Computing_ , vol. 69, no. 2, pp. 117–124, 2009.
  * [41] R. Rabenseifner, “Optimization of collective reduction operations,” 06 2004, pp. 1–9.
  * [42] S. Rajbhandari, J. Rasley, O. Ruwase, and Y. He, “Zero: Memory optimizations toward training trillion parameter models,” 2019.
  * [43] J. Rehr, F. Vila, J. Gardner, L. Svec, and M. Prange, “Scientific computing in the cloud,” _Computing in Science and Engineering_ , vol. 12, pp. 34 – 43, 07 2010.
  * [44] A. Roy, H. Zeng, J. Bagga, G. Porter, and A. C. Snoeren, “Inside the social network’s (datacenter) network,” _Computer Communication Review_ , vol. 45, pp. 123–137, 2015.
  * [45] P. D. Sack, “Scalable collective message-passing algorithms,” Ph.D. dissertation, Champaign, IL, USA, 2011, aAI3503864.
  * [46] P. Sanders, J. Speck, and J. L. Träff, “Two-tree algorithms for full bandwidth broadcast, reduction and scan,” _Parallel Comput._ , vol. 35, no. 12, p. 581–594, Dec. 2009. [Online]. Available: https://doi.org/10.1016/j.parco.2009.09.001
  * [47] R. Thakur, R. Rabenseifner, and W. Gropp, “Optimization of collective communication operations in mpich,” _Int. J. High Perform. Comput. Appl._ , vol. 19, no. 1, p. 49–66, Feb. 2005. [Online]. Available: https://doi.org/10.1177/1094342005051521
  * [48] A. Vishnu, C. Siegel, and J. Daily, “Distributed tensorflow with mpi,” 2016.
  * [49] D. W. Walker and J. J. Dongarra, “Mpi: a standard message passing interface,” _Supercomputer_ , vol. 12, pp. 56–68, 1996.
  * [50] B. Zhang, Y. Ruan, and J. Qiu, “Harp: Collective communication on hadoop,” in _2015 IEEE International Conference on Cloud Engineering_ , 2015, pp. 228–233.