Centralized Network Utility Maximization over Aggregate Flows

TL;DR

This paper presents a decomposition approach for network utility maximization that simplifies complex flow rate problems into aggregate flow problems, with solutions applicable to both single-path and multipath scenarios, using ADMM.

Contribution

It introduces a novel decomposition theorem for NUM problems based on source-destination aggregation, applicable to various utility functions and implementable via ADMM.

Findings

Theorems relating detailed flow rates to aggregate flows for NUM.

A simple apportionment method for $\

Implementation of the decomposition using ADMM.

Abstract

We study a network utility maximization (NUM) decomposition in which the set of flow rates is grouped by source-destination pairs. We develop theorems for both single-path and multipath cases, which relate an arbitrary NUM problem involving all flow rates to a simpler problem involving only the aggregate rates for each source-destination pair. The optimal aggregate flows are then apportioned among the constituent flows of each pair. This apportionment is simple for the case of $\alpha$-fair utility functions. We also show how the decomposition can be implemented with the alternating direction method of multipliers (ADMM) algorithm.

Tables

Table 1. TABLE I: Comparison of ADMM, CP, and Grad. Proj. for Small Graph Example

	ADMM				Gradient Projection				Chambolle-Pock
$N$	$f^{\ast }$	$l_{\max }$	$n_{\text{iter}}$	$t$ (sec)	$f^{\ast }$	$l_{\max }$	$n_{\text{iter}}$	$t$ (sec)	$f^{\ast }$	$l_{\max }$	$n_{\text{iter}}$	$t$ (sec)
10	-92.084	10.000	194	0.0193	-92.085	10.000	431	1.1607	-92.084	10.000	74	0.0030
15	-136.800	10.000	207	0.0259	-136.809	10.000	500	1.2534	-136.800	10.000	112	0.0050
20	-182.002	10.000	304	0.0446	-182.002	10.000	590	1.6405	-182.002	10.000	258	0.0130
25	-243.806	10.000	296	0.0508	-243.806	10.000	1288	3.4499	-243.806	10.000	225	0.0130
30	-289.040	10.000	296	0.0574	-289.040	10.000	2029	5.9341	-289.040	10.000	256	0.0170

Table 2. TABLE II: Parameters for Small Graph Example

	ADMM		Grad Proj	Chambolle-Pock
$N$	$r$	pct	$\alpha$	$\sigma$	$\tau$	$\theta$
10	20	${10}^{-4}$	$1.07\times {10}^{-2}$	1.0	0.020	1.0
15	20	${10}^{-4}$	$1.07\times {10}^{-2}$	1.0	0.015	1.0
20	20	${10}^{-4}$	$1.31\times {10}^{-2}$	1.0	0.015	1.0
25	20	${10}^{-4}$	$6.17\times {10}^{-3}$	1.0	0.015	1.0
30	20	${10}^{-4}$	$6.19\times {10}^{-3}$	1.0	0.013	1.0

Table 3. TABLE III: Comparison of ADMM and CP for Large Graph Example

	ADMM						Chambolle-Pock
$N$	$r$	pct	$f^{\ast }$	$l_{\max }$	$n_{\text{iter}}$	$t$ (sec)	$\sigma$	$\tau$	$\theta$	$f^{\ast }$	$l_{\max }$	$n_{\text{iter}}$	$t$ (sec)
50	40	${10}^{-4}$	-1326.781	10.000	100	0.0358	10.0	$3.0\times {10}^{-4}$	1.0	-1326.780	10.000	70	0.0120
75	40	${10}^{-4}$	-2002.522	10.001	162	0.0777	10.0	$2.0\times {10}^{-4}$	1.0	-2002.522	10.000	314	0.0760
100	40	${10}^{-4}$	-2589.978	10.000	179	0.1093	10.0	$5.0\times {10}^{-4}$	0.1	-2589.978	10.000	745	0.2280
125	40	${10}^{-4}$	-3333.174	10.007	208	0.1529	10.0	$4.9\times {10}^{-4}$	0.1	-3333.174	10.000	843	0.3119

Equations

$$R_{l,i}=\left\{\begin{tabular}[]{ll}$1$,&traffic of class $i$ passes through link $l$\\ $0$,&otherwise.\end{tabular}\right.$$

$$\underset{\{u_{i,k}\}}{\mathrm{maximize}}$$

$$\sum _{i=1}^N\sum _{k=1}^{K_i}{R}_{l,i}{u}_{i,k}\le c_l,l=1,\dots ,L$$

$$\underset{\{x_i\}}{\mathrm{maximize}}$$

$$\sum _{i=1}^N{R}_{l,i}{x}_i\le c_l,l=1,\dots ,L$$

$$(f_1\diamond f_2)(x)=\underset{\{(x_1,x_2):x_1+x_2=x\}}{\sup }{f}_1(x_1)+f_2(x_2).$$

$$f_i(x_i)=\underset{\begin{array}{c}\{u_{i,1},\dots ,u_{i,K_i}\}:\\ {\sum }_k{u}_{i,k}=x_i\end{array}}{\sup }\sum _k{f}_{i,k}(u_{i,k})=(f_{i,1}\diamond \cdots \diamond f_{i,K_i})(x_i)$$

$$\underset{\{u_{i,k}\}}{\mathrm{maximize}}\sum _{k=1}^{K_i}{f}_{i,k}(u_{i,k}),\text{subject to}\sum _{k=1}^{K_i}{u}_{i,k}=x_i^{\ast }$$

$$f_i=f_{i,1}\diamond \cdots \diamond f_{i,K_i}=(\sum _k{f}_{i,k}^{\ast }{)}^{\ast }.$$

$$(f_{i,k}^{\ast }{)}^{\prime }=f_{i,k}^{\prime -1}.$$

$$u_{i,k}^{\ast }=g_{i,k}^{\prime }(f_i^{\prime }(x_i^{\ast })),\forall i,k$$

$${\mathcal{L}}_1(\mathbf{u},\mathbf{\rho })=\sum _i\sum _k{f}_{i,k}(u_{i,k})-\sum _{l=1}^L{\rho }_l(\sum _i\sum _k{R}_{l,i}{u}_{i,k}-c_l).$$

$$\sum _i{R}_{l,i}\sum _k{u}_{i,k}$$

$$\mathbf{\rho }$$

$${\rho }_l(\sum _i{R}_{l,i}\sum _k{u}_{i,k}-c_l)$$

$$u_{i,k}$$

$$h_i=\sum _k{h}_{i,k}=\sum _k{g}_{i,k}^{\prime }=g_i^{\prime }$$

$${\mathcal{L}}_2(\mathbf{x},\mathbf{\lambda })=\sum _i{f}_i(x_i)-\sum _l{\lambda }_l(\sum _i{R}_{l,i}{x}_i-c_l).$$

$$\sum _i{R}_{l,i}{x}_i$$

$$\mathbf{\lambda }$$

$${\lambda }_l(\sum _i{R}_{l,i}{x}_i-c_l)$$

$$x_i$$

$$u_{i,k}=h_{i,k}({\mathbf{\lambda }}^{\ast T}{\mathbf{r}}_i)=h_{i,k}(f_i^{\prime }(x_i^{\ast })),\forall i,k.$$

$$\sum _k{u}_{i,k}=\sum _k{h}_{i,k}({\mathbf{\lambda }}^{\ast T}{\mathbf{r}}_i)=h_i({\mathbf{\lambda }}^{\ast T}{\mathbf{r}}_i)=x_i^{\ast },\forall i$$

$$u_{i,k}^{\ast }=\frac{w_{i,k}}{{\sum }_{k^{\prime }}{w}_{i,k^{\prime }}}{x}_i^{\ast }.$$

$$u_{i,k}=\frac{w_{i,k}}{w_i}{x}_i^{\ast }$$

$$\sum _i{w}_i\frac{{\hat{x}}_i-x_i^{\ast }}{x_i^{\ast }}\le 0.$$

$$0\ge \sum _i\frac{w_i}{x_i^{\ast }}\sum _k({\hat{u}}_{i,k}-u_{i,k})=\sum _i\sum _k{w}_{i,k}\frac{{\hat{u}}_{i,k}-u_{i,k}}{u_{i,k}}$$

$$u_{i,k}^{\ast }=g_{i,k}^{\prime }(f_i^{\prime }(x_i^{\ast }))=\frac{w_{i,k}^{\frac{1}{a+1}}}{{\sum }_{k^{\prime }}{w}_{i,k^{\prime }}^{\frac{1}{a+1}}}{x}_i^{\ast }.$$

$$[{\mathbf{S}}_{i}]_{l,j}=\left\{\begin{tabular}[]{ll}$1$,&\begin{tabular}[]{@{}l@{}}traffic on the $j$th path of class $i$\\ passes through link $l$\end{tabular}\\ $0$,&otherwise.\end{tabular}\right.$$

$$\underset{\{u_{i,j,k}\}}{\mathrm{maximize}}$$

$$\sum _{i=1}^N\sum _{j=1}^J\sum _{k=1}^{K_i}[{\mathbf{S}}_i{]}_{l,j}{u}_{i,j,k}\le c_l,\forall l$$

$$\underset{\{x_{i,j}\}}{\mathrm{maximize}}$$

$$\sum _{i=1}^N\sum _{j=1}^J[{\mathbf{S}}_i{]}_{l,j}{x}_{i,j}\le c_l,\forall l$$

$$\phi_{i,k}({\mathbf{u}}_{i,k})=\left\{\begin{tabular}[]{ll}$f_{i,k}({\mathbf{1}}_{J}^{T}{\mathbf{u}}_{i,k})$,&${\mathbf{u}}_{i,k}\in\mathbb{R}^{J}_{+}$\\ $-\infty$,&otherwise\end{tabular}\right.$$

Full text

I Introduction

The last two decades have seen a great deal of research in network utility maximization (NUM) [1] [2] [3], which has cast light on traditional networking protocols [4] and has facilitated the design of promising future protocols [5] as well. Most NUM researchers have focused on developing distributed solutions to various utility maximization problems. These distributed solutions, which follow nicely from dual decompositions [2], are ideal for internets, in which cooperation among flow sources cannot be assumed, and minimal communication between links and nodes is desired. In recent years there has been growing interest in the software defined networking (SDN) paradigm, in which data and control planes are separated [6] [7]. In this framework, certain network functions such as flow control, congestion control, and throughput optimization may be assigned to a central controller. Central control is feasible for closed networks, such as in data centers [8] or communication satellite networks [9].

In some networks with central control, the number of flows $K$ may be much larger than the number of source-destination pairs $N$. For example, the Iridium satellite network employs 66 satellites and facilitates tens of thousands of flows [10]. A similar phenomenon may occur in small data centers. In this paper, we study a primal decomposition in which the set of flows is grouped into flow classes, each corresponding to a source-destination pair. Many congestion control algorithms inherently group flows by source-destination pair [11] and several related primal decompositions have been studied, for example in [2]. However, because the source-destination decomposition is only applicable to centralized control, it has received little attention. Given the recent popularity of SDN, however, the decomposition may prove to be beneficial. To this end, we develop a comprehensive theory of the source-destination decomposition in this paper. (We also discuss briefly in Section V a potential benefit of this decomposition in a network with “semi-distributed” control.) We derive theorems that decompose a NUM problem with $K$ variables into one with only $N$ variables, followed by an allocation problem which apportions the aggregate rate for each class among the class’s constituent flows. In some cases, this apportionment is simple, In other cases, the alternating direction method of multipliers (ADMM) algorithm can exploit the decompostion numerically.

The remainder of the paper is organized as follows. In Section II we present the aggregate flow decomposition and the main results relating the original NUM problem to the simpler aggregate flow problem. This analysis is extended to the multipath case in Section III. In Section IV we discuss numerical algorithms, which exploit the aggregate flow decomposition. Finally, in Section V, we conclude the paper.

II Optimization Over Aggregate Flows

Consider a communication network with $M$ nodes and $L$ links (edges). Let $N$ be the number of source-destination pairs in use among all flows. Then $N\leq M(M-1)$. Number the source-destination pairs $1,\dots,N$ and call the set of flows in pair $i$ the $i$th flow class. Let $K_{i}$ be the number of flows in class $i$. Then the total number of flows is $K=\sum_{i=1}^{N}K_{i}$. Let $u_{i,k}$ be the rate of the $k$th flow in class $i$. Finally, define the $L\times N$ binary routing matrix ${\mathbf{R}}$ as

[TABLE]

Note that all flows within a class follow the same path. Consider the following utility maximization problem.

[TABLE]

where $f_{i,k}$ is a utility function for the $k$th flow in class $i$, $c_{l}$ is the capacity of link $l$, and the constraints imply that no link is overloaded. Next, let $x_{i}=\sum_{k=1}^{K_{i}}u_{i,k}$ be the aggregate rate of class $i$ and consider the aggregate flow utility maximization problem

[TABLE]

where $f_{i}$ is an aggregate utility function for class $i$. The domains of the utility functions $f_{i,k}$ and aggregate utility functions $f_{i}$ are not stated here, but are usually subsets of $\mathbb{R}_{+}$ as negative flow rates are not allowed.

II-A Decomposition by Supremal Convolutions

Definition 1

Let the functions $f_{1}$ and $f_{2}$ be concave and proper on $\mathbb{R}^{n}$. The supremal convolution of $f_{1}$ and $f_{2}$ is

[TABLE]

The supremal convolution of the concave functions $f_{1}$ and $f_{2}$ is simply the negative infimal convolution of the convex functions $-f_{1}$ and $-f_{2}$. By [12, Theorem 5.4], $f_{1}\diamond f_{2}$ is concave. Observe that problem (2) is equivalent to problem (1) when each $f_{i}$ is the $K_{i}$-fold supremal convolution

[TABLE]

and for each flow class $i$, the optimal subflow rates solve

[TABLE]

where $x_{i}^{*}$ is the solution of problem (2) (provided it exists) with $f_{i}$ defined as above. Thus, when the utilities are concave and proper, problem (1) can be decomposed into an aggregate optimization and $N$ optimizations over the subflows as long as the supremal convolutions can be calculated. This decomposition lends itself to parallel implementations, as the $N$ subproblems are independent.

Let $f^{*}(y)=\inf_{x}(xy-f(x))$ denote the concave Fenchel conjugate. From [12, Theorem 16.4], the conjugate supremal convolution is $(f_{i,1}\diamond\dots\diamond f_{i,K_{i}})^{*}=\sum_{k}f_{i,k}^{*}$. Thus the concave closure of the supremal convolution is $(\sum_{k}f_{i,k}^{*})^{*}$. By [12, Corollary 20.1.1], if the $f_{i,k}$’s are closed and $\cap_{k}\text{relint}(\text{dom}f_{i,k}^{*})\neq\emptyset$, then $f_{i,1}\diamond\dots\diamond f_{i,K_{i}}$ is closed, so

[TABLE]

II-B Decomposition with Functions of Legendre Type

Definition 2

A pair $(f,{\mathcal{D}})$ is of Legendre type if ${\mathcal{D}}$ is a nonempty open convex set, $f$ is a strictly concave differentiable function on ${\mathcal{D}}$ and $\lim_{n\to\infty}\|\nabla f(x_{n})\|=+\infty$ for any sequence $\{x_{n}\}$ in ${\mathcal{D}}$ converging to a boundary point of ${\mathcal{D}}$.

Although the Legendre type property applies to a pair $(f,{\mathcal{D}})$, we will refer to a function $f$ as being of Legendre type when $(f,\text{int}(\text{dom}f))$ is of Legendre type. Note that when $\text{dom}f=\mathbb{R}_{++}$, the last condition in Definition 2 is equivalent to $\lim_{x\downarrow 0}f^{\prime}(x)=+\infty$.

The (concave) Legendre conjugate (or Legendre transform) [12] [13] of a pair $(f,{\mathcal{D}})$, where ${\mathcal{D}}\subset\mathbb{R}$ is open and $f$ is differentiable on ${\mathcal{D}}$, is the pair $(g,{\mathcal{E}})$ where $g(y)=yf^{\prime-1}(y)-f(f^{\prime-1}(y))$ and ${\mathcal{E}}$ is the image of ${\mathcal{D}}$ under $f^{\prime}$.

From [12, Theorem 26.5], if $f_{i,k}$ is closed, ${\mathcal{D}}=\text{int}(\text{dom}f_{i,k})$, and ${\mathcal{D}}^{*}=\text{int}(\text{dom}f_{i,k}^{*})$, then $(f_{i,k},{\mathcal{D}})$ is of Legendre type if and only if $(f^{*}_{i,k},{\mathcal{D}}^{*})$ is of Legendre type. When these pairs are of Legendre type, $(f_{i,k}^{*},{\mathcal{D}}^{*})$ is the Legendre conjugate of $(f_{i,k},{\mathcal{D}})$, which is the Legendre conjugate of $(f_{i,k}^{*},{\mathcal{D}}^{*})$, so conjugation is involutory: $(f^{**}_{i,k},{\mathcal{D}}^{**})=(f_{i,k},{\mathcal{D}})$, and

[TABLE]

Note that if the $f_{i,k}$’s are closed and Legendre type, with domain ${\mathcal{D}}$ and $\cap_{k}\text{relint}(\text{dom}f_{i,k}^{*})\neq\emptyset$, then $\sum_{k}f^{*}_{i,k}$ is Legendre type (and therefore differentiable), and the supremal convolution $(\sum_{k}f^{*}_{i,k})^{*}$ is closed and Legendre type.

Let the $f_{i,k}$’s have domain ${\mathcal{D}}=\mathbb{R}_{++}$. Since there are no equality constraints in problems (1) and (2) and the inequality constraints are all affine, Slater’s condition guarantees strong duality for each problem as long as a feasible point exists in the relative interior of the problem domain [13, Sec. 5.2.3], which is $\mathbb{R}_{++}^{K}$ for problem (1) and $\mathbb{R}_{++}^{N}$ for (2). Clearly, setting all optimization variables to a small $\epsilon>0$ yields such a point, so strong duality holds for both problems. With Legendre type functions, problem (1) is strictly concave with convex feasible region and has a unique solution. Thus there is a unique primal-dual optimal pair satisfying the Karush-Kuhn-Tucker (KKT) conditions for problem (1) with Legendre-type utility functions.

Theorem 1

Let the functions $f_{i,k}$ be closed, concave, and Legendre type with domain $\mathbb{R}_{++}$ and $\cap_{k}\text{relint}(\text{dom}f_{i,k}^{*})\neq\emptyset$. For each $i,k$, let $g_{i,k}=f_{i,k}^{*}$, $g_{i}=\sum_{k}g_{i,k}$, and $f_{i}=g_{i}^{*}$. Let $\{x_{i}^{*}\}$ be a primal solution to problem (2) with this definition of $\{f_{i}\}$. Then (1) has unique primal solution

[TABLE]

*and the corresponding dual solutions of (1) and (2) are equal.

Note that $f_{i}=f_{i,1}\diamond\dots\diamond f_{i,K_{i}}$.*

Proof. First note that $f_{i}$ and $g_{i,k}$ are Legendre type and therefore differentiable. Let $h_{i,k}=g^{\prime}_{i,k}$ for each $i,k$. From (4) we have $h_{i,k}=f^{\prime-1}_{i,k}$. The Lagrangian for problem (1) is

[TABLE]

The KKT sufficient conditions for optimality of (1) are

[TABLE]

where ${\mathbf{r}}_{i}$ is the $i$th column of ${\mathbf{R}}$. Condition (9) is equivalent to $\partial{\mathcal{L}}_{1}/\partial u_{i,k}=0$. Now, let

[TABLE]

for each $i$ and consider problem (2) with $f_{i}=g_{i}^{*}$ Since $f_{i}$ is Legendre type, it is strictly concave and thus $\{x_{i}^{*}\}$ is the unique solution to problem (2). Since, in addition, $f_{i}$ is closed, we can use (4) to get $h_{i}=f_{i}^{\prime-1}$. The Lagrangian is

[TABLE]

The KKT conditions for problem (2) are thus

[TABLE]

Next, let $({\mathbf{x}}^{*},{\boldsymbol{\lambda}}^{*})$ be the primal-dual solution for problem (2) and set ${\boldsymbol{\rho}}={\boldsymbol{\lambda}}^{*}$. Then condition (7) is immediately satisfied. Next let

[TABLE]

Then condition (9) is satisfied and using $h_{i}=f_{i}^{\prime-1}=\sum_{k}h_{i,k}$, and (10), we have

[TABLE]

which ensures $\sum_{k}u_{i,k}^{*}=x_{i}^{*}$, and therefore conditions (6) and (8) are satisfied. Finally, since the image of $h_{i,k}$ is $\mathbb{R}_{++}$, $u_{i,k}>0$ for each $i,k$. ∎

II-C Examples

Here we apply the aggregate flow decomposition to some example utility maximization problems. Theorem 1 can be applied to utility functions belonging to the family of $\alpha$-fair functions [14], while the decomposition using supremal convolutions must be used for more general problems.

II-C1 Weighted Logarithm Utilities

Let $f_{i,k}(u_{i,k})=w_{i,k}\log u_{i,k}$ with $w_{i,k}\geq 0$, so that the overall utility is a sum of weighted logarithms of individual flows. These functions belong to the class of $\alpha$-fair utilities with $\alpha=0$ and are appealing as they yield proportionally fair rate allocations [14]. They are also clearly Legendre type so Theorem 1 can be used. The Legendre conjugates can be calculated using (4). We have $f_{i,k}^{\prime}(u)=w_{i,k}/u$, $g^{\prime}_{i,k}(v)=w_{i,k}/v$, and $g^{\prime}_{i}(v)=\sum_{k}w_{i,k}/v$. Next $f_{i}^{\prime}(x)=g_{i}^{\prime-1}(x)=\sum_{k}w_{i,k}/x$ and $f_{i}(x)=\sum_{k}w_{i,k}\log x$. Finally, using (5)

[TABLE]

Note that all of the above functions have domain $\mathbb{R}_{++}$. From (11), the optimized aggregate flows should be apportioned to the subflows in proportion to their weights. Note that the utility function of problem (2) is also a sum of weighted logarithms, where the $i$th weight is the sum weight of the $i$th class.

The weighted logarithm case can also be proven using proportional fairness [15]. Let $\{x_{i}^{*}\}$ be the solution to problem (2) with $f_{i}(x_{i})=w_{i}\log x_{i}$ and let the subflow rates be

[TABLE]

where $\{w_{i,k}\}$ are any non-negative weights such that $\sum_{k}w_{i,k}=w_{i}$. From [15], the unique solution to (2) $\{x_{i}^{*}\}$ is such that the rates per unit charge are proportionally fair. That is, if $\{\hat{x}_{i}\}$ is any other set of rates then

[TABLE]

Now let $\{\hat{u}_{i,k}\}$ be any set of subflow rates not equal to those found by (12) and let $\hat{x}_{i}=\sum_{k}\hat{u}_{i,k}$. From (12) we have $w_{i}/x_{i}^{*}=w_{i,k}/u_{i,k}$ for all $i,k$. From (13) we have

[TABLE]

So the proportionally allocated solution is such that the rates per unit charge are proportionally fair. Thus it is the unique solution to (1) with $f_{i,k}(u_{i,k})=w_{i,k}\log u_{i,k}$.

II-C2 Weighted Power Utilities (Negative Exponent)

As another example, let $f_{i,k}(u_{i,k})=-w_{i,k}u_{i,k}^{-a}$ with $a\geq 1$. These functions are also part of the $\alpha$-fair family. When $a=1$, the allocation satisfies minimum potential delay fairness and as $a\to+\infty$, the allocation is max-min fair [14]. The utilities are also of Legendre type and we can use Theorem 1. We have $f_{i,k}^{\prime}(u)=aw_{i,k}u^{-(a+1)}$ and $g^{\prime}_{i,k}(v)=(aw_{i,k}/v)^{\frac{1}{a+1}}$. Next we have $g^{\prime}_{i}(y)=\sum_{k}g^{\prime}_{i,k}(y)=f_{i}^{\prime-1}(y)$. Thus $f_{i}^{\prime}(x)=(a/x^{a+1})(\sum_{k}w_{i,k}^{\frac{1}{a+1}})^{a+1}$ and the optimum subflow rates are

[TABLE]

The utility functions for problem (2) are $f_{i}(x_{i})=-x_{i}^{-a}(\sum_{k}w_{i,k}^{\frac{1}{a+1}})^{a+1}$. Again, all of the above functions have domain $\mathbb{R}_{++}$.

II-C3 Quadratic Utilities

Quadratic functions are not of Legendre type and are not necessarily increasing on $\mathbb{R}_{+}$, rendering them unsuitable for use as utility functions. However, the aggregate flow decomposition can be useful when implementing a gradient projection algorithm. In a gradient projection algorithm, steepest ascent iterations are followed by projections onto the feasible set [16]. Such a projection is a quadratic program (QP) that can be simplified by decomposing with supremal convolutions.

Let $\{z_{i,k}\}$ be the set of variables obtained after an iteration of steepest ascent for problem (1). This set must be projected onto the routing polytope $\{\{u_{i,k}\}:\sum_{i}\sum_{k}R_{l,i}u_{i,k}\leq c_{l},l=1,\dots,L\}$. The projection QP is problem (1) with quadratic utility $f_{i,k}(u_{i,k})=-\frac{1}{2}(u_{i,k}-z_{i,k})^{2}$, and domain $\{u_{i,k}\geq 0\}$. (In this section all functions are equal to $-\infty$ outside their domains). Note that $f_{i,k}$ is not Legendre type. However $f_{i,k}$ is closed, concave, and proper on $\mathbb{R}$ and has conjugate $f_{i,k}^{*}(y)=-\frac{1}{2}y^{2}+z_{i,k}y$, with domain $\{y\leq z_{i,k}\}$. Thus any point less than $z_{i,\min}=\min_{k}z_{i,k}$ lies in $\text{relint}(\text{dom}f_{i,k}^{*})$ for all $(i,k)$, and therefore (3) can be used to find $f_{i}$. The conjugate aggregate utility is $f_{i}^{*}(y)=\sum_{k}f_{i,k}^{*}(y)=-K_{i}y^{2}/2+\bar{z}_{i}y$, with domain $\{y\leq z_{i,\min}\}$, where $\bar{z}_{i}=\sum_{k}z_{i,k}$. The aggregate function is obtained by conjugating $f_{i}^{*}$, which yields $f_{i}(x_{i})=-\frac{1}{2K_{i}}(x_{i}-\bar{z}_{i})^{2}$, with domain $\{x_{i}\geq\bar{z}_{i}-K_{i}z_{i,\min}\}$. Finally, for each class $i$, the subflows minimize $\sum_{k}\frac{1}{2}(u_{i,k}-z_{i,k})^{2}$ subject to $\sum_{k}u_{i,k}=x^{*}_{i}$ and $u_{i,k}\geq 0$ for each $k$, where $x_{i}^{*}$ is the solution to the aggregate problem. (Thus $x_{i}^{*}\geq\bar{z}_{i}-K_{i}z_{i,\min}$.) The subflow problem is strictly convex and has unique solution $u_{i,k}^{*}=z_{i,k}+\frac{1}{K_{i}}(x_{i}^{*}-\bar{z}_{i})$.

II-C4 Piecewise Linear Utilities

Piecewise linear functions are important as they are often used as approximations of functions that are difficult to work with analytically or are incompletely known [17]. In this case Theorem 1 is not applicable but supremal convolutions can be calculated using (3). Let $f_{i,k}$ be concave and piecewise linear with non-negative breakpoints $0=c_{1}<c_{2}<\dots<c_{B}$ and corresponding non-negative slopes $m_{1}>m_{2}>\dots>m_{B}=0$, and let $f_{i,k}(c_{1})=f_{i,k}(0)=0$ and $f_{i,k}(x)=-\infty$ for $x<0$. (The number of breakpoints $B$ need not be the same for all utilities.) Then $f_{i,k}$ is closed and from [17, Sec. 8F], the conjugate of $f_{i,k}$ is also concave and piecewise linear with breakpoints $0=m_{B}<m_{B-1}<\dots<m_{1}$ and corresponding slopes $c_{B}>c_{B-1}>\dots>c_{1}$, and $f_{i,k}^{*}(m_{1})=0$. That is, the breakpoints of $f_{i,k}^{*}$ are the slopes of $f_{i,k}$ and the slopes of $f_{i,k}^{*}$ are the breakpoints of $f_{i,k}$. Finally, $\text{dom}f_{i,k}^{*}=\mathbb{R}_{+}$ and thus (3) can be used.

The aggregate utility function $f_{i}$ can be found with the following prescription: For each $f_{i,k}$, find $f_{i,k}^{*}$ by exchanging breakpoints and slopes, as described above. Sum these conjugates to find the conjugate of the aggregate utility $f_{i}^{*}=\sum_{k}f_{i,k}^{*}$. Thus $f_{i}^{*}$ is piecewise linear and concave as well. Finally, exchange slopes and breakpoints of $f_{i}^{*}$ to arrive at $f_{i}$.

Therefore, the piecewise-linear problem, which is a linear program (LP) in $K=\sum_{i}K_{i}$ variables, can be decomposed into one LP in $N$ variables, followed by $N$ parallel sub-LP’s, the $i$th sub-LP having $K_{i}$ variables.

III Extension to Multipath Case

Now suppose that for each flow class $i$, traffic may be split into subflows and routed over multiple paths. (Here, a subflow refers to that portion of a flow routed over a certain path, as opposed to a constituent flow of a flow class). Let $J$ be the number of paths and assume $J$ is the same for all flow classes. For each class $i$, define the $L\times J$ routing matrix ${\mathbf{S}}_{i}$ as

[TABLE]

and let the overall $L\times NJ$ routing matrix be ${\mathbf{R}}=[{\mathbf{S}}_{1},\dots,{\mathbf{S}}_{N}]$. Finally, let $u_{i,j,k}$ be the rate on the $j$th path of flow $k$ of class $i$. The multipath utility maximization problem is

[TABLE]

Letting $x_{i,j}=\sum_{k}u_{i,j,k}$ be the aggregate rate on path $j$ of class $i$, the aggregate flow problem in the multipath case is

[TABLE]

For these problems, explicit constraints for non-negativity of the throughputs are added because, for example $\sum_{j}u_{i,j,k}$ can be non-negative even with some negative subflows. Unlike the single-path case, neither problem is strictly convex.

III-A Multipath Supremal Convolution Decomposition

Let the $f_{i,k}$’s be concave and proper on $\mathbb{R}$ with $\text{dom}f_{i,k}\subset\mathbb{R}_{+}$. For each $i,k$ pair, define

[TABLE]

where ${\mathbf{u}}_{i,k}=[u_{i,1,k},\dots,u_{i,J,k}]^{T}$. Then $\phi_{i,k}$ is concave and proper on $\mathbb{R}^{J}$ (but not strictly concave, even if $f_{i,k}$ is) and problem (14) is equivalent to

[TABLE]

subject to the link load constraints of (14). Define the aggregate flow problem by

[TABLE]

with the link load constraints of (15). Here $\phi_{i}({\mathbf{x}}_{i})$ is a function from $\mathbb{R}^{J}$ to $\mathbb{R}$ and ${\mathbf{x}}_{i}=[x_{i,1},\dots,x_{i,J}]^{T}$. Similar to the argument in Section II-A, problem (17) is equivalent to problem (16) with concave aggregate functions $\phi_{i}=\phi_{i,1}\diamond\dots\diamond\phi_{i,K_{i}}$ if for each flow class $i$, the optimal subflow rates solve the problem

[TABLE]

where ${\mathbf{x}}_{i}^{*}$ is the solution of problem (17) (provided it exists) with $\phi_{i}$ defined as above. If the $f_{i,k}$’s are closed and $\cap_{k}\text{relint}(\text{dom}f_{i,k}^{*})\neq\emptyset$, then it can be shown that the same is true of the $\phi_{i,k}$’s and the aggregate functions can be found using $\phi_{i}=(\sum_{k}\phi_{i,k}^{*})^{*}$ where $\phi^{*}({\mathbf{y}})=\inf_{{\mathbf{u}}}({\mathbf{u}}^{T}{\mathbf{y}}-\phi({\mathbf{u}}))$.

III-B Multipath Legendre-Type Case

As in the single-path case, when $\text{dom}f_{i,k}=\mathbb{R}_{++}$, strong duality of (14) and (15) follows from Slater’s condition. However, in the multipath case, neither problem is strictly convex and uniqueness of the solutions cannot be guaranteed.

Theorem 2

Let $f_{i,k}$, $g_{i,k}$, $f_{i}$, and $g_{i}$ satisfy the conditions of Theorem 1 for all $i,k$. Let $\{x_{i,j}^{*}\}$ be a primal solution to (15) with this definition of $\{f_{i}\}$. Then a solution of the following constrained system of linear equations

[TABLE]

(where $\bar{x}_{i}^{*}=\sum_{j}x_{i,j}^{*}$), is a primal solution to (14). Furthermore, if ${\boldsymbol{\lambda}}^{*}$ is the dual solution to (15) corresponding to the link load constraints and $\mu_{i,j}^{*}$ is the dual solution to (15) corresponding to the non-negativitiy constraint of $x_{i,j}$, then the dual solution to (14) corresponding to the link load constraints is ${\boldsymbol{\rho}}^{*}={\boldsymbol{\lambda}}^{*}$ and the dual solution to (14) corresponding to the non-negativity constraint of $u_{i,j,k}$ is $\sigma_{i,j,k}^{*}=\mu_{i,j}^{*}$, for each $k$.

Proof. Let $h_{i,k}=g^{\prime}_{i,k}=f^{\prime-1}_{i,k}$ for each $i,k$ and let $\bar{u}_{i,k}=\sum_{j}u_{i,j,k}$. The Lagrangian for problem (14) is

[TABLE]

Setting the derivative with respect to $u_{i,j,k}$ to zero gives

[TABLE]

Thus, the following seven relations constitute the KKT conditions for problem (14):

[TABLE]

Next, turning to problem (15) with $f_{i}=g_{i}^{*}$, let $h_{i}=\sum_{k}h_{i,k}=\sum_{k}g^{\prime}_{i,k}=g_{i}^{\prime}$ for each $i$ and let $\bar{x}_{i}=\sum_{j}x_{i,j}$. The Lagrangian is

[TABLE]

The KKT conditions for problem (15) are

[TABLE]

Now set ${\boldsymbol{\rho}}={\boldsymbol{\lambda}}^{*}$ and $\sigma_{i,j,k}=\mu_{i,j}^{*}$ for all $i,j,k$ and let $\{u_{i,j,k}^{*}\}$ be a solution to (18). Then it can be seen that all KKT conditions (20)–(26) are satisfied, and furthermore, using (19)

[TABLE]

which ensures $\sum_{j,k}u_{i,j,k}^{*}=\sum_{j}x_{i,j}^{*}$. ∎

Note that when $J=1$, the second equation of problem (18) reduces to (5) which guarantees that the first equation and the non-negativity condition hold.

The subflow allocation problem given by (18) can be decomposed into $N$ parallel problems (one for each class $i$). Let ${\mathbf{u}}_{i,k}=[u_{i,1,k},\dots,u_{i,J,k}]^{T}$, ${\mathbf{u}}_{i}=[{\mathbf{u}}_{i,1}^{T},\dots,{\mathbf{u}}_{i,K_{i}}^{T}]^{T}$, ${\mathbf{x}}_{i}=[x_{i,1},\dots,x_{i,J}]^{T}$, and define the matrices ${\mathbf{A}}_{i}={\mathbf{1}}_{K_{i}}^{T}\otimes{\mathbf{I}}_{J}$ and ${\mathbf{B}}_{i}={\mathbf{I}}_{K_{i}}\otimes{\mathbf{1}}_{J}^{T}$. Then the $i$th optimal subflow vector ${\mathbf{u}}_{i}^{*}$ solves

[TABLE]

where ${\mathbf{g}}_{i}=[g_{i,1}^{\prime}(f_{i}^{\prime}(\bar{x}_{i}^{*})),\dots,g_{i,K_{i}}^{\prime}(f_{i}^{\prime}(\bar{x}_{i}^{*}))]^{T}$.

Note that the only component of (31) that depends on the utility functions is ${\mathbf{g}}_{i}$. As an example, for the case of weighted logarithm utilities (see Section II-C1) with $f_{i,k}(u)=w_{i,k}\log u$, we have ${\mathbf{g}}_{i}=(\bar{x}_{i}^{*}/\bar{w}_{i}){\mathbf{w}}_{i}$, where ${\mathbf{w}}_{i}=[w_{i,1},\dots,w_{i,K_{i}}]^{T}$, and $\bar{w}_{i}=\sum_{k}w_{i,k}$.

IV Aggregate Decomposition with ADMM

Here we show that the alternating direction method of multipliers (ADMM) algorithm [18] can inherently decompose problem (1) into an optimization over aggregate flows and $N$ parallel optimizations over the constituent flows. Assume that the utility functions $f_{i,k}$ have domain $\mathbb{R}_{+}$.

IV-A ADMM Algorithm

To apply ADMM to problem (1) we recast it as

[TABLE]

where ${\mathbf{x}}=[x_{1},\dots,x_{N}]^{T}$ and ${\mathbf{y}}=[y_{1},\dots,y_{L}]^{T}$. The function $h$ indicates that the links are not overloaded. That is $h({\mathbf{y}})=0$ if ${\mathbf{y}}\leq{\mathbf{c}}$, and $+\infty$ otherwise, with ${\mathbf{c}}=[c_{1},\dots,c_{L}]^{T}$. The augmented Lagrangian for problem (32) is

[TABLE]

where ${\mathbf{s}}=\left[\sum_{k}u_{1,k},\dots,\sum_{k}u_{N,k}\right]^{T}$. Here ${\boldsymbol{\lambda}}\in\mathbb{R}^{N}$ and ${\boldsymbol{\rho}}\in\mathbb{R}^{L}$ are the dual variables and $r$ is the penalty parameter. The ADMM method involves repeated minimizations of ${\mathcal{L}}_{r}$ over $({\mathbf{u}},{\mathbf{y}})$, and then ${\mathbf{x}}$. The minimizer with respect to ${\mathbf{x}}$ is ${\mathbf{x}}={\mathbf{A}}^{-1}{\mathbf{b}}$ where ${\mathbf{A}}={\mathbf{I}}+{\mathbf{R}}^{T}{\mathbf{R}}$ and ${\mathbf{b}}={\mathbf{s}}+{\mathbf{R}}^{T}{\mathbf{y}}+({\boldsymbol{\lambda}}+{\mathbf{R}}^{T}{\boldsymbol{\rho}})/r$. Minimization of ${\mathcal{L}}_{r}$ with respect to ${\mathbf{y}}$ is decoupled from that of ${\mathbf{u}}$ and is achieved by a simple projection of ${\mathbf{R}}{\mathbf{x}}-{\boldsymbol{\rho}}/r$ onto the box $\{{\mathbf{y}}:{\mathbf{y}}\leq{\mathbf{c}}\}$. Finally, minimization with respect to ${\mathbf{u}}$ involves $N$ parallel minimizations of the form

[TABLE]

IV-B Numerical Examples

Here we apply the ADMM algorithm with the aggregate flow decomposition to a few example cases and compare performance against a gradient projection algorithm and the primal-dual algorithm of Chambolle and Pock [19]. We examine two sample graphs. The first, which we call the small graph, shown in Figure 1, has $M=6$ nodes and $L=14$ links. The maximum number of source-destination pairs is $N_{\max}=30$. The second example graph, the large graph is shown in Figure 2. This represents the topology of the Iridium low earth orbit satellite constellation [10], and includes $M=66$ satellites (nodes) and $L=192$ links (the actual topology changes as satellites enter and exit polar regions). The dark nodes in Figure 2 represent satellites linked to ground stations which connect flows to terrestrial networks. Thus, we assume all flows either originate or terminate at one of these nodes. The resulting maximum number of source-destination pairs is $N_{\max}=750$. For both example graphs, we set all link capacities to 10 units. Thus ${\mathbf{c}}=10\cdot{\mathbf{1}}_{L}$. In all examples, the number of flows in any flow class (source-destination pair) is uniformly distributed between 10 and 20. The total number of flows is thus $15N$ on average. Finally, the path (route) for each source-destination pair is found using Dijkstra’s algorithm.

IV-B1 ADMM with Weighted Logarithm Utilities

We let the utility function for the $k$th flow of class $i$ be $f_{i,k}(u_{i,k})=w_{i,k}\log u_{i,k}$ with weights $w_{i,k}$ chosen uniformly from $(0,1)$ and solve the optimization problem with ADMM. The minimizer of the augmented Lagrangian with respect to the individual flows $\{u_{i,k}\}$ is found by solving (33) for each $i$ which gives

[TABLE]

where $\psi_{i}=\lambda_{i}-rx_{i}$ and $\bar{w}_{i}=\sum_{k}w_{i,k}$. The ADMM iteration is then given by

[TABLE]

where $[\cdot]^{+}$ represents projection onto the box $\{{\mathbf{y}}:{\mathbf{y}}\leq{\mathbf{c}}\}$, and ${\mathbf{s}}^{(n+1)}=[\sum_{k}u_{1,k}^{(n+1)},\dots,\sum_{k}u_{N,k}^{(n+1)}]^{T}$.

IV-B2 Gradient Projection Optimizer

We compare the ADMM algorithm with a simple gradient projection optimizer. The gradient projection optimizer utilizes Theorem 1 with aggregate utilities $f_{i}(x_{i})=\bar{w}_{i}\log x_{i}$ and optimal subflow rates $u_{i,k}^{*}=w_{i,k}x_{i}^{*}/\bar{w}_{i}$ (see Section II-C1). The update rule for the aggregate problem is

[TABLE]

where $\alpha>0$ is a step size and $P_{{\mathbf{R}}}$ is the function which projects onto the routing polytope $\{{\mathbf{x}}:{\mathbf{R}}{\mathbf{x}}\leq{\mathbf{c}},{\mathbf{x}}\geq\bf 0\}$. In all examples that follow, $P_{{\mathbf{R}}}$ , which solves a QP, is implemented using the CVXOPT QP solver [20].

IV-B3 Chambolle-Pock Optimizer

Problem (1) can be solved with the Chambolle-Pock (CP) algorithm by writing it as

[TABLE]

where ${\mathbf{u}}$ is the concatenation of the $N$ subflow rate vectors $\{{\mathbf{u}}_{i}\}$ and $g$ is the indicator function of the box $\{{\mathbf{y}}:{\mathbf{y}}\leq{\mathbf{c}}\}$. The $L\times K$ matrix ${\mathbf{Q}}$ is defined by

[TABLE]

The algorithm requires evaluation of the proximal operators [12] of $\sigma g^{*}$ and $\tau f$ where $f({\mathbf{u}})=-\sum_{i,k}f_{i,k}(u_{i,k})$, $g^{*}$ is the convex conjugate of $g$, and $\sigma$ and $\tau$ are positive constants. Using Moreau’s theorem [12, Theorem 31.5] we get ${\mbox{\rm prox}}_{\sigma g^{*}}({\mathbf{z}})={\mathbf{z}}-\sigma[{\mathbf{z}}/\sigma]^{+}$ (again $[\cdot]^{+}$ signifies projection onto $\{{\mathbf{y}}:{\mathbf{y}}\leq{\mathbf{c}}\}$). The proximal operator of $\tau f$ is

[TABLE]

The algorithm consists of the following iteration

[TABLE]

where $\theta\in[0,1]$.

IV-B4 Algorithm Comparison with Small Graph

Each iteration of ADMM contains three sparse (0-1 matrix)-vector multiplies with ${\mathbf{R}}$ and ${\mathbf{R}}^{T}$ and one $N\times N$ set of linear equations with the same coefficient matrix ${\mathbf{A}}$. The CP iterations contain two multiplications with the sparse 0-1 matrices ${\mathbf{Q}}$ and ${\mathbf{Q}}^{T}$. Finally, each iteration of the gradient projection algorithm solves a QP (with sparse constraint matrix ${\mathbf{G}}=[{\mathbf{R}}^{T},-{\mathbf{I}}_{N}]^{T}$). Thus, the gradient algorithm has the highest per-iteration cost, followed by ADMM and CP.

The algorithms’ performances are summarized in Table I for various numbers of source-destination pairs $N$. For each algorithm, the converged objective value $f^{*}$ is shown, along with the maximum link load $l_{\max}$, the number of iterations $n_{\text{iter}}$, and the optimization time $t$. In Table II the algorithm parameters are listed, including the ADMM penalty parameter $r$, the gradient projection step-size $\alpha$, and the percent threshold (pct). This value is used as the stopping criterion for ADMM (i.e., when the augmented Lagrangian changes by less than pct percent, stop). Also shown are the three CP parameters $\sigma$, $\tau$, and $\theta$. The gradient projection step-sizes and CP parameters are individually tuned for fastest convergence. The optimization times are averaged over 10 runs (with identical random number generator seeds). All simulations were performed using Python/Numpy, and the projection step in the gradient projection algorithm uses the CVXOPT QP solver (which in turn uses the CHOLMOD sparse Cholesky solver). The optimization time of ADMM and CP is plotted versus $N$ in Figure 3. From Table I, the number of CP iterations required for this example is consistently less than the number of ADMM iterations. As it has a lower per-iteration cost, the convergence time of CP is lower. The gradient projection algorithm has the highest per-iteration cost as well as the largest number of iterations, and thus converges slowest. Note that, without the aggregate flow decomposition (Theorem 1), the gradient projection optimizer would be far slower.

IV-B5 Algorithm Comparison with Large Graph

Next, we repeat the experiment using the large graph. Table III shows the results along with the algorithm parameters. The gradient projection algorithm has been omitted as its convergence times are far greater than ADMM and CP. The optimization times for the ADMM and CP algorithms are plotted in Figure 4. In this example, as $N$ increases, the number of ADMM iterations grows slower than the number of CP iterations. Thus, although the CP per-iteration cost is lower, the larger number of iterations for large $N$ renders CP slower than ADMM.

V Conclusion

We have shown that for many types of utilities, the solution to a $K$-flow NUM problem can be found by solving a simpler $N$-variable problem. This principle holds for both single-path and multipath NUM problems. The results of this paper have applicability for software-defined networks in which a controller must solve the global NUM problem. These results can also be beneficial for networks consisting of several hub-spoke clusters. For example, with $N$ clusters and $K_{i}$ sources in the $i$th cluster, the problem (2) can be substituted for problem (1). This simpler problem could then be solved in a distributed manner by the hub nodes, which would in turn allocate subflow rates to the spoke nodes.

VI Acknowledgements

The first author wishes to thank Dr. Joseph Yadegar from UtopiaCompression Corporation and Dr. You Lu of Google for insightful conversations regarding this work. This research was partly supported by the United States Air Force under contract number FA9453-14-C-0060. The views and conclusions contained herein are those of the authors and should not be interpreted as necessarily representing the official policies or endorsements, either expressed or implied, of the United States Air Force.