Lyapunov optimization

Last updated October 22, 2022

This article describes Lyapunov optimization for dynamical systems. It gives an example application to optimal control in queueing networks.

Introduction

Lyapunov optimization refers to the use of a Lyapunov function to optimally control a dynamical system. Lyapunov functions are used extensively in control theory to ensure different forms of system stability. The state of a system at a particular time is often described by a multi-dimensional vector. A Lyapunov function is a nonnegative scalar measure of this multi-dimensional state. Typically, the function is defined to grow large when the system moves towards undesirable states. System stability is achieved by taking control actions that make the Lyapunov function drift in the negative direction towards zero.

Lyapunov drift is central to the study of optimal control in queueing networks. A typical goal is to stabilize all network queues while optimizing some performance objective, such as minimizing average energy or maximizing average throughput. Minimizing the drift of a quadratic Lyapunov function leads to the backpressure routing algorithm for network stability, also called the max-weight algorithm.^[1]^[2] Adding a weighted penalty term to the Lyapunov drift and minimizing the sum leads to the drift-plus-penalty algorithm for joint network stability and penalty minimization.^[3]^[4]^[5] The drift-plus-penalty procedure can also be used to compute solutions to convex programs and linear programs.^[6]

Lyapunov drift for queueing networks

Consider a queueing network that evolves in discrete time with normalized time slots $t\in \{0,1,2,\ldots \}.$ Suppose there are $N$ queues in the network, and define the vector of queue backlogs at time $t$ by:

Q(t)=(Q_{1}(t),\ldots ,Q_{N}(t))

Quadratic Lyapunov functions

For each slot $t,$ define:

L(t)={\frac {1}{2}}\sum _{i=1}^{N}Q_{i}(t)^{2}

This function is a scalar measure of the total queue backlog in the network. It is called quadratic Lyapunov function on the queue state. Define the Lyapunov drift as the change in this function from one slot to the next:

\Delta L(t)=L(t+1)-L(t)

Bounding the Lyapunov drift

Suppose the queue backlogs change over time according to the following equation:

Q_{i}(t+1)=\max \left\{Q_{i}(t)+a_{i}(t)-b_{i}(t),0\right\}

where $a_{i}(t)$ and $b_{i}(t)$ are arrivals and service opportunities, respectively, in queue $i$ on slot $t.$ This equation can be used to compute a bound on the Lyapunov drift for any slot t:

Q_{i}(t+1)^{2}=\left(\max \left\{Q_{i}(t)+a_{i}(t)-b_{i}(t),0\right\}\right)^{2}\leqslant \left(Q_{i}(t)+a_{i}(t)-b_{i}(t)\right)^{2}

Rearranging this inequality, summing over all $i,$ and dividing by 2 leads to:

\Delta L(t)\leqslant B(t)+\sum _{i=1}^{N}Q_{i}(t)(a_{i}(t)-b_{i}(t))\qquad (Eq.1)

where:

B(t)={\frac {1}{2}}\sum _{i=1}^{N}\left(a_{i}(t)-b_{i}(t)\right)^{2}

Suppose the second moments of arrivals and service in each queue are bounded, so that there is a finite constant $B>0$ such that for all $t$ and all possible queue vectors $Q(t)$ the following property holds:

\mathbb {E} [B(t)|Q(t)]\leqslant B

Taking conditional expectations of (Eq. 1) leads to the following bound on the conditional expected Lyapunov drift:

\mathbb {E} [\Delta L(t)|Q(t)]\leqslant B+\sum _{i=1}^{N}Q_{i}(t)\mathbb {E} [a_{i}(t)-b_{i}(t)|Q(t)]\qquad (Eq.2)

A basic Lyapunov drift theorem

In many cases, the network can be controlled so that the difference between arrivals and service at each queue satisfies the following property for some real number $\varepsilon >0$ :

\mathbb {E} [a_{i}(t)-b_{i}(t)|Q(t)]\leqslant -\varepsilon

If the above holds for the same epsilon for all queues $i,$ all slots $t,$ and all possible vectors $Q(t),$ then (Eq. 2) reduces to the drift condition used in the following Lyapunov drift theorem. The theorem below can be viewed as a variation on Foster's theorem for Markov chains. However, it does not require a Markov chain structure.

Theorem (Lyapunov Drift).^[5]^[7] Suppose there are constants

B\geqslant 0,\varepsilon >0

such that for all

t

and all possible vectors

Q(t)

the conditional Lyapunov drift satisfies:

\mathbb {E} [\Delta L(t)|Q(t)]\leqslant B-\varepsilon \sum _{i=1}^{N}Q_{i}(t).

Then for all slots

t>0

the time average queue size in the network satisfies:

{\frac {1}{t}}\sum _{\tau =0}^{t-1}\sum _{i=1}^{N}\mathbb {E} [Q_{i}(\tau )]\leqslant {\frac {B}{\varepsilon }}+{\frac {\mathbb {E} [L(0)]}{\varepsilon t}}.

Proof. Taking expectations of both sides of the drift inequality and using the law of iterated expectations yields:

\mathbb {E} [\Delta L(t)]\leqslant B-\varepsilon \sum _{i=1}^{N}\mathbb {E} [Q_{i}(t)]

Summing the above expression over $\tau \in \{0,1,\ldots ,t-1\}$ and using the law of telescoping sums gives:

\mathbb {E} [L(t)]-\mathbb {E} [L(0)]\leqslant Bt-\varepsilon \sum _{\tau =0}^{t-1}\sum _{i=1}^{N}\mathbb {E} [Q_{i}(\tau )]

Using the fact that $L(t)$ is non-negative and rearranging the terms in the above expression proves the result.

Lyapunov optimization for queueing networks

Consider the same queueing network as in the above section. Now define $p(t)$ as a network penalty incurred on slot $t.$ Suppose the goal is to stabilize the queueing network while minimizing the time average of $p(t).$ For example, to stabilize the network while minimizing time average power, $p(t)$ can be defined as the total power incurred by the network on slot t.^[8] To treat problems of maximizing the time average of some desirable reward $r(t),$ the penalty can be defined $p(t)=-r(t).$ This is useful for maximizing network throughput utility subject to stability.^[3]

To stabilize the network while minimizing the time average of the penalty $p(t),$ network algorithms can be designed to make control actions that greedily minimize a bound on the following drift-plus-penalty expression on each slot $t$ :^[5]

\Delta L(t)+Vp(t)

where $V$ is a non-negative weight that is chosen as desired to affect a performance tradeoff. A key feature of this approach is that it typically does not require knowledge of the probabilities of the random network events (such as random job arrivals or channel realizations). Choosing $V=0$ reduces to minimizing a bound on the drift every slot and, for routing in multi-hop queueing networks, reduces to the backpressure routing algorithm developed by Tassiulas and Ephremides.^[1]^[2] Using $V>0$ and defining $p(t)$ as the network power use on slot $t$ leads to the drift-plus-penalty algorithm for minimizing average power subject to network stability developed by Neely.^[8] Using $V>0$ and using $p(t)$ as the negative of an admission control utility metric leads to the drift-plus-penalty algorithm for joint flow control and network routing developed by Neely, Modiano, and Li.^[3]

A generalization of the Lyapunov drift theorem of the previous section is important in this context. For simplicity of exposition, assume $p(t)$ is bounded from below:

p(t)\geqslant p_{\min }\quad \forall t\in \{0,1,2,...\}

For example, the above is satisfied with $p_{\min }=0$ in cases when the penalty $p(t)$ is always non-negative. Let $p^{*}$ represent a desired target for the time average of $p(t).$ Let $V$ be a parameter used to weight the importance of meeting the target. The following theorem shows that if a drift-plus-penalty condition is met, then the time average penalty is at most O(1/V) above the desired target, while average queue size is O(V). The $V$ parameter can be tuned to make time average penalty as close to (or below) the target as desired, with a corresponding queue size tradeoff.

Theorem (Lyapunov Optimization). Suppose there are constants

\varepsilon >0,V,B\geqslant 0,

and

p^{*}

such that for all

t

and all possible vectors

Q(t)

the following drift-plus-penalty condition holds:

\mathbb {E} [\Delta L(t)+Vp(t)|Q(t)]\leqslant B+Vp^{*}-\varepsilon \sum _{i=1}^{N}Q_{i}(t)

Then for all

t>0

the time average penalty and time average queue sizes satisfy:

{\frac {1}{t}}\sum _{\tau =0}^{t-1}\mathbb {E} [p(\tau )]\leqslant p^{*}+{\frac {B}{V}}+{\frac {\mathbb {E} [L(0)]}{Vt}}

{\frac {1}{t}}\sum _{\tau =0}^{t-1}\sum _{i=1}^{N}\mathbb {E} [Q_{i}(\tau )]\leqslant {\frac {B+V(p^{*}-p_{\min })}{\varepsilon }}+{\frac {\mathbb {E} [L(0)]}{\varepsilon t}}

Proof. Taking expectations of both sides of the posited drift-plus-penalty and using the law of iterated expectations we have:

\mathbb {E} [\Delta L(t)]+V\mathbb {E} [p(t)]\leqslant B+Vp^{*}-\varepsilon \sum _{i=1}^{N}\mathbb {E} [Q_{i}(t)]

Summing the above over the first $t$ slots and using the law of telescoping sums gives:

{\begin{aligned}\mathbb {E} [L(t)]-\mathbb {E} [L(0)]+V\sum _{\tau =0}^{t-1}\mathbb {E} [p(\tau )]&\leqslant (B+Vp^{*})t-\varepsilon \sum _{\tau =0}^{t-1}\sum _{i=1}^{N}\mathbb {E} [Q_{i}(\tau )]\\-\mathbb {E} [L(0)]+V\sum _{\tau =0}^{t-1}\mathbb {E} [p(\tau )]&\leqslant (B+Vp^{*})t&&{\text{Since }}L(t),Q_{i}(t)\geqslant 0\\V\sum _{\tau =0}^{t-1}\mathbb {E} [p(\tau )]&\leqslant p^{*}Vt+Bt+\mathbb {E} [L(0)]\end{aligned}}

Dividing by $Vt$ and rearranging terms proves the time average penalty bound. A similar argument proves the time average queue size bound.

Related Research Articles

The Beer–Lambert law, also known as Beer's law, the Lambert–Beer law, or the Beer–Lambert–Bouguer law relates the attenuation of light to the properties of the material through which the light is travelling. The law is commonly applied to chemical analysis measurements and used in understanding attenuation in physical optics, for photons, neutrons, or rarefied gases. In mathematical physics, this law arises as a solution of the BGK equation.

In differential geometry, a subject of mathematics, a symplectic manifold is a smooth manifold, $, equipped with a closed nondegenerate differential 2-form, called the symplectic form. The study of symplectic manifolds is called symplectic geometry or symplectic topology. Symplectic manifolds arise naturally in abstract formulations of classical mechanics and analytical mechanics as the cotangent bundles of manifolds. For example, in the Hamiltonian formulation of classical mechanics, which provides one of the major motivations for the field, the set of all possible configurations of a system is modeled as a manifold, and this manifold's cotangent bundle describes the phase space of the system.$

<span class="mw-page-title-main">Fokker–Planck equation</span> Partial differential equation

In statistical mechanics, the Fokker–Planck equation is a partial differential equation that describes the time evolution of the probability density function of the velocity of a particle under the influence of drag forces and random forces, as in Brownian motion. The equation can be generalized to other observables as well.

In condensed matter physics, Bloch's theorem states that solutions to the Schrödinger equation in a periodic potential take the form of a plane wave modulated by a periodic function. The theorem is named after the physicist Felix Bloch, who discovered the theorem in 1929. Mathematically, they are written

Various types of stability may be discussed for the solutions of differential equations or difference equations describing dynamical systems. The most important type is that concerning the stability of solutions near to a point of equilibrium. This may be discussed by the theory of Aleksandr Lyapunov. In simple terms, if the solutions that start out near an equilibrium point $stay near forever, then is Lyapunov stable . More strongly, if is Lyapunov stable and all solutions that start out near converge to, then is asymptotically stable. The notion of exponential stability guarantees a minimal rate of decay, i.e., an estimate of how quickly the solutions converge. The idea of Lyapunov stability can be extended to infinite-dimensional manifolds, where it is known as structural stability, which concerns the behavior of different but "nearby" solutions to differential equations. Input-to-state stability (ISS) applies Lyapunov notions to systems with inputs.$

Transmittance of the surface of a material is its effectiveness in transmitting radiant energy. It is the fraction of incident electromagnetic power that is transmitted through a sample, in contrast to the transmission coefficient, which is the ratio of the transmitted to incident electric field.

In mathematics, the Poisson summation formula is an equation that relates the Fourier series coefficients of the periodic summation of a function to values of the function's continuous Fourier transform. Consequently, the periodic summation of a function is completely defined by discrete samples of the original function's Fourier transform. And conversely, the periodic summation of a function's Fourier transform is completely defined by discrete samples of the original function. The Poisson summation formula was discovered by Siméon Denis Poisson and is sometimes called Poisson resummation.

The Drude model of electrical conduction was proposed in 1900 by Paul Drude to explain the transport properties of electrons in materials. Basically, Ohm's law was well established and stated that the current J and voltage V driving the current are related to the resistance R of the material. The inverse of the resistance is known as the conductance. When we consider a metal of unit length and unit cross sectional area, the conductance is known as the conductivity, which is the inverse of resistivity. The Drude model attempts to explain the resistivity of a conductor in terms of the scattering of electrons by the relatively immobile ions in the metal that act like obstructions to the flow of electrons.

In physics, the Hamilton–Jacobi equation, named after William Rowan Hamilton and Carl Gustav Jacob Jacobi, is an alternative formulation of classical mechanics, equivalent to other formulations such as Newton's laws of motion, Lagrangian mechanics and Hamiltonian mechanics. The Hamilton–Jacobi equation is particularly useful in identifying conserved quantities for mechanical systems, which may be possible even when the mechanical problem itself cannot be solved completely.

In mathematics, more specifically in the field of analytic number theory, a Landau–Siegel zero or simply Siegel zero, named after Edmund Landau and Carl Ludwig Siegel, is a type of potential counterexample to the generalized Riemann hypothesis, on the zeros of Dirichlet L-functions associated to quadratic number fields. Roughly speaking, these are possible zeros very near to $s = 1$ .

In mathematics, a Kloosterman sum is a particular kind of exponential sum. They are named for the Dutch mathematician Hendrik Kloosterman, who introduced them in 1926 when he adapted the Hardy–Littlewood circle method to tackle a problem involving positive definite diagonal quadratic forms in four as opposed to five or more variables, which he had dealt with in his dissertation in 1924.

In chaos theory, the correlation integral is the mean probability that the states at two different times are close:

<span class="mw-page-title-main">Vertex model</span>

A vertex model is a type of statistical mechanics model in which the Boltzmann weights are associated with a vertex in the model. This contrasts with a nearest-neighbour model, such as the Ising model, in which the energy, and thus the Boltzmann weight of a statistical microstate is attributed to the bonds connecting two neighbouring particles. The energy associated with a vertex in the lattice of particles is thus dependent on the state of the bonds which connect it to adjacent vertices. It turns out that every solution of the Yang–Baxter equation with spectral parameters in a tensor product of vector spaces $yields an exactly-solvable vertex model.$

Anatoly Alexeyevich Karatsuba was a Russian mathematician working in the field of analytic number theory, p-adic numbers and Dirichlet series.

Defect types include atom vacancies, adatoms, steps, and kinks that occur most frequently at surfaces due to the finite material size causing crystal discontinuity. What all types of defects have in common, whether surface or bulk defects, is that they produce dangling bonds that have specific electron energy levels different from those of the bulk. This difference occurs because these states cannot be described with periodic Bloch waves due to the change in electron potential energy caused by the missing ion cores just outside the surface. Hence, these are localized states that require separate solutions to the Schrödinger equation so that electron energies can be properly described. The break in periodicity results in a decrease in conductivity due to defect scattering.

In mathematics, the Levi-Civita field, named after Tullio Levi-Civita, is a non-Archimedean ordered field; i.e., a system of numbers containing infinite and infinitesimal quantities. Each member $can be constructed as a formal series of the form$

In queueing theory, a discipline within the mathematical theory of probability, the backpressure routing algorithm is a method for directing traffic around a queueing network that achieves maximum network throughput, which is established using concepts of Lyapunov drift. Backpressure routing considers the situation where each job can visit multiple service nodes in the network. It is an extension of max-weight scheduling where each job visits only a single service node.

In thermal quantum field theory, the Matsubara frequency summation is the summation over discrete imaginary frequencies. It takes the following form

In the mathematical theory of probability, the drift-plus-penalty method is used for optimization of queueing networks and other stochastic systems.

In mathematics, the injective tensor product of two topological vector spaces (TVSs) was introduced by Alexander Grothendieck and was used by him to define nuclear spaces. An injective tensor product is in general not necessarily complete, so its completion is called the completed injective tensor products. Injective tensor products have applications outside of nuclear spaces. In particular, as described below, up to TVS-isomorphism, many TVSs that are defined for real or complex valued functions, for instance, the Schwartz space or the space of continuously differentiable functions, can be immediately extended to functions valued in a Hausdorff locally convex TVS $with out any need to extend definitions from real/complex-valued functions to -valued functions.$

References

1 2 L. Tassiulas and A. Ephremides, "Stability Properties of Constrained Queueing Systems and Scheduling Policies for Maximum Throughput in Multihop Radio Networks, IEEE Transactions on Automatic Control, vol. 37, no. 12, pp. 1936-1948, Dec. 1992.
1 2 L. Tassiulas and A. Ephremides, "Dynamic Server Allocation to Parallel Queues with Randomly Varying Connectivity," IEEE Transactions on Information Theory, vol. 39, no. 2, pp. 466-478, March 1993.
1 2 3 M. J. Neely, E. Modiano, and C. Li, "Fairness and Optimal Stochastic Control for Heterogeneous Networks," Proc. IEEE INFOCOM, March 2005.
↑ L. Georgiadis, M. J. Neely, and L. Tassiulas, "Resource Allocation and Cross-Layer Control in Wireless Networks," Foundations and Trends in Networking, vol. 1, no. 1, pp. 1-149, 2006.
1 2 3 M. J. Neely. Stochastic Network Optimization with Application to Communication and Queueing Systems, Morgan & Claypool, 2010.
↑ M. J. Neely, "Distributed and Secure Computation of Convex Programs over a Network of Connected Processors," DCDIS Conf, Guelph, Ontario, July 2005
↑ E. Leonardi, M. Mellia, F. Neri, and M. Ajmone Marsan, "Bounds on Average Delays and Queue Size Averages and Variances in Input-Queued Cell-Based Switches", Proc. IEEE INFOCOM, 2001.
1 2 M. J. Neely, "Energy Optimal Control for Time Varying Wireless Networks," IEEE Transactions on Information Theory, vol. 52, no. 7, pp. 2915-2934, July 2006.

Primary Sources

M. J. Neely. Stochastic Network Optimization with Application to Communication and Queueing Systems, Morgan & Claypool, 2010.

This page is based on this Wikipedia article
Text is available under the CC BY-SA 4.0 license; additional terms may apply.
Images, videos and audio are available under their respective licenses.

[tass-radio-nets-1] 1 2 L. Tassiulas and A. Ephremides, "Stability Properties of Constrained Queueing Systems and Scheduling Policies for Maximum Throughput in Multihop Radio Networks, IEEE Transactions on Automatic Control, vol. 37, no. 12, pp. 1936-1948, Dec. 1992.

[tass-server-allocation-2] 1 2 L. Tassiulas and A. Ephremides, "Dynamic Server Allocation to Parallel Queues with Randomly Varying Connectivity," IEEE Transactions on Information Theory, vol. 39, no. 2, pp. 466-478, March 1993.

[neely-fairness-infocom05-3] 1 2 3 M. J. Neely, E. Modiano, and C. Li, "Fairness and Optimal Stochastic Control for Heterogeneous Networks," Proc. IEEE INFOCOM, March 2005.

[now-4] L. Georgiadis, M. J. Neely, and L. Tassiulas, "Resource Allocation and Cross-Layer Control in Wireless Networks," Foundations and Trends in Networking, vol. 1, no. 1, pp. 1-149, 2006.

[sno-text-5] 1 2 3 M. J. Neely. Stochastic Network Optimization with Application to Communication and Queueing Systems, Morgan & Claypool, 2010.

[neely-dcdis-6] M. J. Neely, "Distributed and Secure Computation of Convex Programs over a Network of Connected Processors," DCDIS Conf, Guelph, Ontario, July 2005

[leonardi-7] E. Leonardi, M. Mellia, F. Neri, and M. Ajmone Marsan, "Bounds on Average Delays and Queue Size Averages and Variances in Input-Queued Cell-Based Switches", Proc. IEEE INFOCOM, 2001.

[neely-energy-it-8] 1 2 M. J. Neely, "Energy Optimal Control for Time Varying Wireless Networks," IEEE Transactions on Information Theory, vol. 52, no. 7, pp. 2915-2934, July 2006.

[1]

[2]

[3]

[4]

[5]

[6]

[7]

[8]