Communication lower bounds and optimal algorithms for programs that reference arrays -- Part 1

Introduction

Algorithms have two costs: computation (e.g., arithmetic) and communication, i.e., moving data between levels of a memory hierarchy or between processors over a network. Communication costs (measured in time or energy per operation) already greatly exceed computation costs, and the gap is growing over time following technological trends [FOSC, ComputerPerformanceNRC]. Thus it is important to design algorithms that minimize communication, and if possible attain communication lower bounds. In this work, we measure communication cost in terms of the number of words moved (a bandwidth cost), and will not discuss other factors like per-message latency, congestion, or costs associated with noncontiguous data. Our goal here is to establish new lower bounds on the communication cost of a much broader class of algorithms than possible before, and when possible describe how to attain these lower bounds.

Communication lower bounds have been a subject of research for a long time. Hong and Kung [hongkung] used an approach called pebbling to establish lower bounds for $\Theta(N^{3})$ matrix multiplication and other algorithms. Irony, Tiskin and Toledo [ITT04] proved the result for $\Theta(N^{3})$ matrix multiplication in a different, geometric way, and extended the results both to the parallel case and the case of using redundant copies of the data. In [BallardDemmelHoltzSchwartz11] this geometric approach was further generalized to include any algorithm that “geometrically resembled” matrix multiplication in a sense to be made clear later, but included most direct linear algebra algorithms (dense or sparse, sequential or parallel), and some graph algorithms as well. Of course lower bounds alone are not algorithms, so a great deal of additional work has gone into developing algorithms that attain these lower bounds, resulting in many faster algorithms as well as remaining open problems (we discuss attainability in Section LABEL:sec:attain).

Our geometric approach, following [ITT04, BallardDemmelHoltzSchwartz11], works as follows. We have a set $\mathcal{Z}$ of arithmetic operations to perform, and the amount of data available locally, i.e., without any communication, is $M$ words. For example, $M$ could be the cache size. Suppose we can upper bound the number of (useful) arithmetic operations that we can perform with just this data; call this bound $F$. Letting $|\cdot|$ denote the cardinality of a set, if the total number of arithmetic operations that we need to perform is $|\mathcal{Z}|$ (e.g., $|\mathcal{Z}|=N^{3}$ multiply-adds in the case of dense matrix multiplication on one processor), then we need to refill the cache at least $|\mathcal{Z}|/F$ times in order to perform all the operations. Since refilling the cache has a communication cost of moving $O(M)$ words (e.g., writing at most $M$ words from cache back to slow memory, and reading at most $M$ new words into cache from slow memory), the total communication cost is $\Omega(M\cdot|\mathcal{Z}|/F)$ words moved. This argument is formalized in Section 4.

The most challenging part of this argument is determining $F$. Our approach, described in more detail in Sections 2–3, builds on the work in [ITT04, BallardDemmelHoltzSchwartz11]; in those papers, the algorithm is modeled geometrically using the iteration space of the loop nest, as sketched in the following example.

Our approach is based on a major generalization of [LW49] in [BCCT10] that lets us geometrically model a much larger class of algorithms with an arbitrary number of loops and array expressions involving affine functions of indices.

formalizes the argument that an upper bound on $F$ yields a lower bound on communication of the form $\Omega(M\cdot|\mathcal{Z}|/F)$. Some technical assumptions are required for this to work. Using the same approach as in [BallardDemmelHoltzSchwartz11], we need to eliminate the possibility that an algorithm could do an unbounded amount of work on a fixed amount of data without requiring any communication. We also describe how the bound applies to communication in a variety of computer architectures (sequential, parallel, heterogeneous, etc.).

summarizes our results, and outlines the contents of Part 2 of this paper. Part 2 will discuss how to compute lower bounds more efficiently and will include more cases where optimal algorithms are possible, including a discussion of loop dependencies.

This completes the outline of the paper. We note that it is possible to omit the detailed proofs in Sections 3 and 5 on a first reading; the rest of the paper is self-contained.

To conclude this introduction, we apply our theory to two examples (revisited later), and show how to derive communication-optimal sequential algorithms.

yielding the well-known blocked algorithm where the innermost three loops multiply $M^{1/2}$–by–$M^{1/2}$ blocks of $A$ and $B$ and update a block of $C$. We note that to have all three blocks fit in fast memory simultaneously, they would have to be slightly smaller than $M^{1/2}$–by–$M^{1/2}$ by a constant factor. We will address this constant factor and others, which are important in practice, in Part 2 of this work (see Section LABEL:sec:conclusions). $\triangle$

Other examples appearing later include matrix-vector multiplication, tensor contractions, the direct $N$–body algorithm, database join, and computing matrix powers $A^{k}$.

Geometric Model

We begin by reviewing the geometric model of matrix multiplication introduced in [ITT04], describe how it was extended to more general linear algebra algorithms in [BallardDemmelHoltzSchwartz11], and finally show how to generalize it to the class of programs considered in this work.

We want a bound $F\geq|E|$, where $E\subseteq\mathcal{Z}$ is any set of lattice points $\mathcal{I}$ representing operations that can be performed just using data available in fast memory of size $M$. Let $E_{A}$, $E_{B}$ and $E_{C}$ be projections of $E$ onto faces of the cube representing all the required operands $A(i_{1},i_{3})$, $B(i_{3},i_{2})$ and $C(i_{1},i_{2})$, resp., needed to perform the operations represented by $E$. A special case of [LW49, Theorem 2] gives us the desired bound on $E$:

Finally, since the entries represented by $E_{A},\allowbreak E_{B},\allowbreak E_{C}$ fit in fast memory by assumption (i.e., $|E_{A}|,\allowbreak|E_{B}|,\allowbreak|E_{C}|\leq M$), this yields the desired bound $F\geq|E|$:

Irony et al. [ITT04] applied this approach to obtain a lower bound for matrix multiplication. Ballard et al. [BallardDemmelHoltzSchwartz11] extended this approach to programs of the form

and picked the iteration space $\mathcal{Z}$, arrays $A,B,C$, and binary operations $+_{i_{1},i_{2}}$ and $\ast_{i_{1},i_{2},i_{3}}$ to represent many (dense or sparse) linear algebra algorithms. To make this generalization, Ballard et al. made several observations, which also apply to our case, below.

An instance of the geometric model is an abstract representation of an algorithm, taking the form

(We will be more concrete about what it means to ‘access a variable’ when we introduce our execution model in Section 4.)

Assuming that all variables are accessed in each iteration rules out certain optimizations, or requires us not to count certain inner loop iterations in our lower bound. For example, consider Boolean matrix multiplication, where the inner loop body is $C(i_{1},i_{2})=C(i_{1},i_{2})\vee(A(i_{1},i_{3})\wedge B(i_{3},i_{2}))$. As long as $A(i_{1},i_{3})\wedge B(i_{3},i_{2})$ is false, $C(i_{1},i_{2})$ does not need to be accessed. So we need either to exclude such optimizations, or not count such loop iterations. And once $C(i_{1},i_{2})$ is true, its value cannot change, and no more values of $A(i_{1},i_{3})$ and $B(i_{3},i_{2})$ need to be accessed; if this optimization is performed, then these loop iterations would not occur and not be counted in the lower bound. In general, if there are data-dependent branches in the inner loop (e.g., ‘flip a coin…’), we may not know to which subset $\mathcal{Z}$ of all the loop iterations our model applies until the code has executed. We refer to a later example (database join) and [BallardDemmelHoltzSchwartz11, Section 3.2.1] for more examples and discussion of this assumption.

Following the approach for matrix multiplication, we need to bound $|E|$ where $E$ is a nonempty finite subset of $\mathcal{Z}$. To execute $E$, we need array entries $A_{1}(\phi_{1}(E)),\ldots,\allowbreak A_{m}(\phi_{m}(E))$. Analogous to the analysis of matrix multiplication, we need to bound $|E|$ in terms of the number of array entries needed: $|\phi_{j}(E)|$ entries of $A_{j}$, for $j=\{{1,\ldots,m}\}$.

Suppose there were nonnegative constants $s_{1},\ldots,s_{m}$ such that for any such $E$,

i.e., a generalization of Theorem 2.1. Then just as with matrix multiplication, if the number of entries of each array were bounded by $|\phi_{j}(E)|\leq M$, we could bound

The next section shows how to determine the constants $s_{1},\ldots,s_{m}$ such that inequality (2.3) holds. To give some rough intuition for the result, consider again the special case of Theorem 2.1, where $|E|$ can be interpreted as a volume, and each $|\phi_{j}(E)|$ as an area. Thus one can think of $|E|$ as having units of, say, ${\rm meters}^{3}$, and $|\phi_{j}(E)|$ as having units of ${\rm meters}^{2}$. Thus, for $\prod_{j=1}^{3}|\phi_{j}(E)|^{s_{j}}$ to have units of ${\rm meters}^{3}$, we would expect the $s_{j}$ to satisfy $3=2s_{1}+2s_{2}+2s_{3}$. But if $E$ were a lower dimensional subset, lying say in a plane, then $|E|$ would have units ${\rm meters}^{2}$ and we would get a different constraint on the $s_{j}$. This intuition is formalized in Theorem 3.2 below.

Upper bounds from HBL theory

As discussed above, we aim to establish an upper bound of the form (2.3) for $|E|$, for any finite set $E$ of loop iterations. In this section we formulate and prove Theorem 3.2, which states that such a bound is valid if and only if the exponents $s_{1},\ldots,s_{m}$ satisfy a certain finite set of linear constraints, specified in (3.1) below.

The following theorem provides bounds which are fundamental to our application.

Conversely, if (3.2) holds for $s\in^{m}$, then $s$ satisfies (3.1).

Note that if $s\in[0,\infty)^{m}$ satisfies (3.1) or (3.2), then any $t\in[0,\infty)^{m}$ where each $t_{j}\geq s_{j}$ does too. This is obvious for (3.1); in the case of (3.2), since $E$ is nonempty and finite, then for each $j$, $1\leq|\phi_{j}(E)|<\infty$, so $|\phi_{j}(E)|^{s_{j}}\leq|\phi_{j}(E)|^{t_{j}}$.

The following result demonstrates there is no loss of generality restricting $s\in^{m}$, rather than in $[0,\infty)^{m}$, in Theorem 3.2.

Whenever (3.1) or (3.2) holds for some $s\in[0,\infty)^{m}$, it holds for $t\in^{m}$ where $t_{j}=\min(s_{j},1)$ for $j\in\{{1,\ldots,m}\}$.

2 Generalizations

We will use the terminology HBL datum to refer to either of two types of structures; the first is defined as follows:

An Abelian group HBL datum is a $3$–tuple

where $G$ and each $G_{j}$ are finitely generated Abelian groups, $G$ is torsion-free, each $\phi_{j}\colon G\to G_{j}$ is a group homomorphism, and the notation on the left-hand side always implies that $j\in\{{1,2,\ldots,m}\}$.

Consider an Abelian group HBL datum $({G,({G_{j}}),({\phi_{j}})})$ and $s\in^{m}$. Suppose that

Conversely, if (3.5) holds for $s\in^{m}$, then $s$ satisfies (3.8).

[BCCT10, Theorem 2.4] treats the more general situation in which $G$ is not required to be torsion-free, and establishes the conclusion (3.4) in the weaker form

where the constant $C<\infty$ depends only on $G$ and $\{{\phi_{1},\ldots,\phi_{m}}\}$. The constant $C=1$ established here is optimal. Indeed, consider any single $x\in G$, and for each $j$ define $f_{j}$ to be $1_{\phi_{j}(x)}$, the indicator function for $\phi_{j}(x)$. Then both sides of the inequalities in (3.4) are equal to $1$. In Remark 3.11, we will show that the constant $C$ is bounded above by the number of torsion elements of $G$.

Necessity of (3.3) follows from an argument given in [BCCT10, Theorem 2.4].

On the other hand, for $j\in\{{1,\ldots,m}\}$,

where $A_{j}$ is a finite constant which depends on $\phi_{j}$, on the structure of $H_{j}$, and on the choice of $\{{e_{i}}\}$, but not on $N$. Indeed, it follows from the definition of rank that for each $j$ it is possible to permute the indices $i$ so that for each $i>\operatorname{rank}({\phi_{j}(H)})$ there exist integers $k_{i}$ and $\kappa_{i,l}$ such that

The upper bound (3.7) follows from these relations.

where $A<\infty$ is independent of $N$. By letting $N$ tend to infinity, we conclude that $\operatorname{rank}({H})\leq{\sum_{j=1}^{m}s_{j}\operatorname{rank}({\phi_{j}(H)})}$, as was to be shown. ∎

We show sufficiency in Theorem 3.6 in its full generality is a consequence of the special case in which all of the groups $G_{j}$ are torsion-free.

We are assuming validity of Theorem 3.6 in the torsion-free case. Its conclusion asserts that

Combining these inequalities gives the conclusion of Theorem 3.6. ∎

Our second generalization of Theorem 3.2 is as follows.

Consider a vector space HBL datum $({V,({V_{j}}),({\phi_{j}})})$ and $s\in^{m}$. Suppose that

Thus the conclusion (3.4) of Theorem 3.6 is satisfied. ∎

Consider $({G,({G_{j}}),({\phi_{j}})})$, an HBL datum with torsion, and $s\in^{m}$. If (3.3) holds, then (3.6) holds with $C=|T(G)|$. In particular,

Conversely, if (3.11) holds for $s\in^{m}$, then $s$ satisfies (3.3).

The first inequality is an application of Theorem 3.6. Summation with respect to $t\in T(G)$ gives the required bound.

The factor $|T(G)|$ cannot be improved if the groups $G_{j}$ are torsion free, or more generally if $T(G)$ is contained in the intersection of the kernels of all the homomorphisms $\phi_{j}$; this is seen by considering $E=T(G)$.

3 The polytope 𝒫𝒫\mathcal{P}

The constraints (3.3) and (3.8) are encoded by convex polytopes in $^{m}$.

For any Abelian group HBL datum $({G,({G_{j}}),({\phi_{j}})})$, we denote the set of all $s\in^{m}$ which satisfy (3.3) by $\mathcal{P}({G,({G_{j}}),({\phi_{j}})})$. For any vector space HBL datum $({V,({V_{j}}),({\phi_{j}})})$, we denote the set of all $s\in^{m}$ which satisfy (3.8) by $\mathcal{P}({V,({V_{j}}),({\phi_{j}})})$.

In the course of proving Theorem 3.6 above, we established the following result.

Now we prove Proposition 3.3, i.e., that there was no loss of generality in assuming each $s_{j}\leq 1$.

Proposition 3.3 concerns the case of Theorem 3.2, but in the course of its proof we will show that this result also applies to Theorems 3.6 and 3.10.

We first show the result concerning (3.1), by considering instead the set of inequalities (3.8). Suppose that a vector space HBL datum $({V,({V_{j}}),({\phi_{j}})})$ and $s\in[0,\infty)^{m}$ satisfy (3.8), and suppose that $s_{k}\geq 1$ for some $k$. Define $t\in[0,\infty)^{m}$ by $t_{j}=s_{j}$ for $j\neq k$, and $t_{k}=1$. Pick any subspace $W\leq V$; let $W^{\prime}\coloneqq W\cap\operatorname{ker}({\phi_{k}})$ and let $U$ be a supplement of $W^{\prime}$ in $W$, i.e., $W=U+W^{\prime}$ and $\operatorname{dim}({W})=\operatorname{dim}({U})+\operatorname{dim}({W^{\prime}})$. Since $\operatorname{dim}({\phi_{k}(U)})=\operatorname{dim}({U})$,

since $s$ satisfies (3.8) and $\operatorname{dim}({\phi_{k}(W^{\prime})})=0$, by (3.8) applied to $W^{\prime}$,

Next, we show the result concerning (3.2). Consider the following more general situation. Let $X,X_{1},\ldots,X_{m}$ be sets and $\phi_{j}\colon X\to X_{j}$ be functions for $j\in\{{1,\ldots,m}\}$. Let $s\in[0,\infty)^{m}$ with some $s_{k}\geq 1$, and suppose that $|E|\leq\prod_{j=1}^{m}|\phi_{j}(E)|^{s_{j}}$ for any finite nonempty subset $E\subseteq X$. Fix one such set $E$. For each $y\in\phi_{k}(E)$, let $E_{y}=\phi_{k}^{-1}(y)\cap E$, the preimage of $y$ under \phi_{k}\big{|}_{E}; thus $|\phi_{k}(E_{y})|=1$. By assumption, $|E_{y}|\leq\prod_{j=1}^{m}|\phi_{j}(E_{y})|^{s_{j}}$, so it follows that

Since $E$ can be written as the union of disjoint sets $\bigcup_{y\in\phi_{k}(E)}E_{y}$, we obtain

We claim that this result can be generalized to $\eqref{BL}$ and $\eqref{BLfield}$, based on a comment in [BCCT10, Section 8] that it generalizes in the weaker case of (3.6). ∎

To prove Theorem 3.10, we show in Section 3.4 that if (3.9) holds at each extreme point $s$ of $\mathcal{P}$, then it holds for all $s\in\mathcal{P}$. Then in Section 3.5, we show that when $s$ is an extreme point of $\mathcal{P}$, the hypothesis (3.8) can be restated in a special form. Finally in Section 3.7 we prove Theorem 3.10 (with restated hypothesis) when $s$ is any extreme point of $\mathcal{P}$, thus proving the theorem for all $s\in\mathcal{P}$.

4 Interpolation between extreme points of 𝒫𝒫\mathcal{P}

One multilinear extension of the Riesz-Thörin theorem states the following (see, e.g., [bennettsharpley]).

Suppose that $p_{0}=({p_{j,0}}),p_{1}=({p_{j,1}})\in[1,\infty]^{m}$. Suppose that there exist $A_{0},A_{1}\in[0,\infty)$ such that

For each $\theta\in(0,1)$ define exponents $p_{j,\theta}$ by

Here $\|f_{j}\|_{p}=\|f_{j}\|_{L^{p}(X_{j},\mathcal{A}_{j},\mu_{j})}$.

In the context of Theorem 3.10 with vector space HBL datum $({V,({V_{j}}),({\phi_{j}})})$, we consider the multilinear map

representing the left-hand side in (3.9).

If (3.9) holds for every extreme point of $\mathcal{P}$, then it holds for every $s\in\mathcal{P}$.

5 Critical subspaces and extreme points

Assume a fixed vector space HBL datum $({V,({V_{j}}),({\phi_{j}})})$, and let $\mathcal{P}$ continue to denote the set of all $s\in^{m}$ which satisfy (3.8).

Consider any $s\in^{m}$. A subspace $W\leq V$ satisfying $\operatorname{dim}({W})=\sum_{j=1}^{m}s_{j}\operatorname{dim}({\phi_{j}(W)})$ is said to be a critical subspace; one satisfying $\operatorname{dim}({W})\leq\sum_{j=1}^{m}s_{j}\operatorname{dim}({\phi_{j}(W)})$ is said to be subcritical; and a subspace satisfying $\operatorname{dim}({W})>\sum_{j=1}^{m}s_{j}\operatorname{dim}({\phi_{j}(W)})$ is said to be supercritical. $W$ is said to be strictly subcritical if $\operatorname{dim}({W})<\sum_{j=1}^{m}s_{j}\operatorname{dim}({\phi_{j}(W)})$.

In this language, the conditions (3.8) assert that that every subspace $W$ of $V$, including $\{{0}\}$ and $V$ itself, is subcritical; equivalently, there are no supercritical subspaces. When more than one $m$–tuple $s$ is under discussion, we sometimes say that $W$ is critical, subcritical, supercritical or strictly subcritical with respect to $s$.

The goal of Section 3.5 is to establish the following:

Let $s$ be an extreme point of $\mathcal{P}$. Then some subspace $\{{0}\}<W<V$ is critical with respect to $s$, or $s\in\{{0,1}\}^{m}$.

Note that these two possibilities are not mutually exclusive.

If $s$ is an extreme point of $\mathcal{P}$, and if $i$ is an index for which $s_{i}\notin\{{0,1}\}$, then $\operatorname{dim}({\phi_{i}(V)})\neq 0$.

Suppose $\operatorname{dim}({\phi_{i}(V)})=0$. If $t\in^{m}$ satisfies $t_{j}=s_{j}$ for all $j\neq i$, then $\sum_{j=1}^{m}t_{j}\operatorname{dim}({\phi_{j}(W)})=\sum_{j=1}^{m}s_{j}\operatorname{dim}({\phi_{j}(W)})$ for all subspaces $W\leq V$, so $t\in\mathcal{P}$ as well. If $s_{i}\notin\{{0,1}\}$, then this contradicts the assumption that $s$ is an extreme point of $\mathcal{P}$. ∎

Let $s$ be an extreme point of $\mathcal{P}$. Suppose that no subspace $\{{0}\}<W\leq V$ is critical with respect to $s$. Then $s\in\{{0,1}\}^{m}$.

Suppose to the contrary that for some index $i$, $s_{i}\notin\{{0,1}\}$. If $t\in^{m}$ satisfies $t_{j}=s_{j}$ for all $j\neq i$ and if $t_{i}$ is sufficiently close to $s_{i}$, then $t\in\mathcal{P}$. This again contradicts the assumption that $s$ is an extreme point. ∎

Let $s$ be an extreme point of $\mathcal{P}$. Suppose that there exists no subspace $\{{0}\}<W<V$ which is critical with respect to $s$. Then there exists at most one index $i$ for which $s_{i}\notin\{{0,1}\}$.

Whenever $|\varepsilon|$ is sufficiently small, $t\in^{m}$. Moreover, $V$ remains subcritical with respect to $t$. If $|\varepsilon|$ is sufficiently small, then every subspace $\{{0}\}<W<V$ remains strictly subcritical with respect to $t$, because the set of all parameters $({\operatorname{dim}({W}),\operatorname{dim}({\phi_{1}(W)}),\ldots,\operatorname{dim}({\phi_{m}(W)})})$ which arise, is finite. Thus $t\in\mathcal{P}$ for all sufficiently small $|\varepsilon|$. Therefore $s$ is not an extreme point of $\mathcal{P}$. ∎

Let $s\in^{m}$. Suppose that $V$ is critical with respect to $s$. Suppose that there exists exactly one index $i\in\{{1,2,\ldots,m}\}$ for which $s_{i}\notin\{{0,1}\}$. Then $V$ has a subspace which is supercritical with respect to $s$.

By Lemma 3.19, $\operatorname{dim}({\phi_{i}(V)})>0$. Let $K$ be the set of all indices $k$ for which $s_{k}=1$. The hypothesis that $V$ is critical means that

Since $s_{i}>0$ and $\operatorname{dim}({\phi_{i}(V)})>0$,

Consider the subspace $W\leq V$ defined by

this intersection is interpreted to be $W=V$ if the index set $K$ is empty. $W$ necessarily has positive dimension. Indeed, $W$ is the kernel of the map $\psi\colon V\to\bigoplus_{k\in K}\phi_{k}(V)$, defined by $\psi(x)=({\phi_{k}(x):k\in K})$, where $\bigoplus$ denotes the direct sum of vector spaces. The image of $\psi$ is isomorphic to some subspace of $\bigoplus_{k\in K}\phi_{k}(V)$, a vector space whose dimension $\sum_{k\in K}\operatorname{dim}({\phi_{k}(V)})$ is strictly less than $\operatorname{dim}({V})$. Therefore $\operatorname{ker}({\psi})=W$ has dimension greater than or equal to $\operatorname{dim}({V})-\sum_{k\in K}\operatorname{dim}({\phi_{k}(V)})>0$. Since $\phi_{k}(W)=\{{0}\}$ for all $k\in K$,

Since $s_{i}<1$ and $\operatorname{dim}({W})>0$, $s_{i}\operatorname{dim}({\phi_{i}(W)})$ is strictly less than $\operatorname{dim}({W})$, whence $W$ is supercritical. ∎

Suppose that there exists no critical subspace $\{{0}\}<W<V$. By Lemma 3.20, either $s_{j}\in\{{0,1}\}^{m}$ — in which case the proof is complete — or $V$ is critical. By Lemma 3.21, there can be at most one index $i$ for which $s_{i}\notin\{{0,1}\}$. By Lemma 3.22, for critical $V$, the existence of one single such index $i$ implies the presence of some supercritical subspace, contradicting the main hypothesis of Proposition 3.18. Thus again, $s\in\{{0,1}\}^{m}$. ∎

6 Factorization of HBL data

Let $({V,({V_{j}}),({\phi_{j}})})$ be a vector space HBL datum. To any subspace $W\leq V$ can be associated two HBL data:

Given the vector space HBL datum $({V,({V_{j}}),({\phi_{j}})})$, for any subspace $W\leq V$,

Consider any subspace $U\leq V$ and some $s\in^{m}$ such that both s\in\mathcal{P}({W,({V_{j}}),({\phi_{j}\big{|}_{W}})}) and $s\in\mathcal{P}({V/W,({V_{j}/\phi_{j}(W)}),({[\phi_{j}]})})$. Then

The last inequality is a consequence of the inclusions $\phi_{j}(U\cap W)\subseteq\phi_{j}(U)\cap\phi_{j}(W)$. The last equality is the relation $\operatorname{dim}({A})+\operatorname{dim}({B})=\operatorname{dim}({A+B})+\operatorname{dim}({A\cap B})$, which holds for any subspaces $A,B$ of a vector space. Thus $U$ is subcritical. ∎

Given the vector space HBL datum $({V,({V_{j}}),({\phi_{j}})})$, let $s\in\mathcal{P}({V,({V_{j}}),({\phi_{j}})})$. Let $W\leq V$ be a subspace which is critical with respect to $s$. Then

With Lemma 3.24 in hand, it remains to show that $\mathcal{P}({V,({V_{j}}),({\phi_{j}})})$ is contained in the intersection of the other two polytopes.

Any subspace $U\leq W$ is also a subspace of $V$. $U$ is subcritical with respect to $s$ when regarded as a subspace of $W$, if and only if $U$ is subcritical when regarded as a subspace of $V$. So s\in\mathcal{P}({W,({V_{j}}),({\phi_{j}\big{|}_{W}})}).

Now consider any subspace of $U/W\leq W/V$; we have $W\leq U\leq V$ and $\operatorname{dim}({U/W})=\operatorname{dim}({U})-\operatorname{dim}({W})$. Moreover,

Therefore since $\operatorname{dim}({W})=\sum_{j=1}^{m}s_{j}\operatorname{dim}({\phi_{j}(W)})$,

by the subcriticality of $U$, which holds because $s\in\mathcal{P}({V,({V_{j}}),({\phi_{j}})})$. Thus any $U/W\leq V/W$ is subcritical with respect to $s$, so $s\in\mathcal{P}({V/W,({V_{j}/\phi_{j}(W)}),({[\phi_{j}]})})$ as well. ∎

7 Proof of Theorem 3.10

Recall we are given the vector space HBL datum $({V,({V_{j}}),({\phi_{j}})})$; we prove Theorem 3.10 by induction on the dimension of the ambient vector space $V$. If $\operatorname{dim}({V})=0$ then $V$ has a single element, and the result (3.9) is trivial.

To establish the inductive step, consider any extreme point $s$ of $\mathcal{P}$. According to Proposition 3.18, there are two cases which must be analyzed. We begin with the case in which there exists a critical subspace $\{{0}\}<W<V$, which we prove in the following lemma. We assume that Theorem 3.10 holds for all HBL data for which the ambient vector space has strictly smaller dimension than is the case for the given datum.

Let $({V,({V_{j}}),({\phi_{j}})})$ be a vector space HBL datum, and let $s\in\mathcal{P}(V,({V_{j}}),({\phi_{j}}))$. Suppose that subspace $\{{0}\}<W<V$ is critical with respect to $s$. Then (3.9) holds for this $s$.

Consider any inequality in (3.9). We may assume that none of the exponents equal zero. For if $s_{k}=0$, then $f_{k}(\phi_{k}(x))\leq\|f_{k}\|_{1/s_{k}}$ for all $x$, and therefore

If $\|f_{k}\|_{1/s_{k}}=0$, then (3.9) holds with both sides . Otherwise we divide by $\|f_{k}\|_{1/s_{k}}$ to conclude that $s\in\mathcal{P}({V,({V_{j}}),({\phi_{j}})})$ if and only if $({s_{j}})_{j\neq k}$ belongs to the polytope associated to the HBL datum $({V,({V_{j}})_{j\neq k},({\phi_{j}})_{j\neq k}})$. Thus the index $k$ can be eliminated. This reduction can be repeated to remove all indices which equal zero.

Let $W_{j}\coloneqq\phi_{j}(W)$. By Lemma 3.25, s\in\mathcal{P}({W,({W_{j}}),({\phi_{j}\big{|}_{W}})}). Therefore by the inductive hypothesis, one of the inequalities in (3.9) is

Define $F_{j}\colon V_{j}/W_{j}\to[0,\infty)$ to be the function

This quantity is a function of the coset $x+W_{j}$ alone, rather than of $x$ itself, because for any $z\in W_{j}$,

by virtue of the substitution $y+z\mapsto y$. Moreover,

To prove this, choose one element $x\in V_{j}$ for each coset $x+W_{j}\in V_{j}/W_{j}$. Denoting by $X$ the set of all these representatives,

because the map $X\times W_{j}\ni(x,y)\mapsto x+y\in V_{j}$ is a bijection.

The inductive bound (3.14) can be equivalently written in the more general form

for any $y\in V$, by applying (3.14) to $({\hat{f}_{j}})$ where $\hat{f}_{j}(z)=f_{j}(z+\phi_{j}(y))$.

Denote by $Y\subseteq V$ a set of representatives of the cosets $y+W\in V/W$, and identify $V/W$ with $Y$. Then

This is a set of inequalities of exactly the form (3.9), with $({V,({V_{j}}),({\phi_{j}})})$ replaced by $({V/W,({V_{j}/W_{j}}),({[\phi_{j}]})})$. By Lemma 3.25, $s\in\mathcal{P}({V/W,({V_{j}/W_{j}}),({[\phi_{j}]})})$, and since $\operatorname{dim}({V/W})<\operatorname{dim}({V})$, we conclude directly from the inductive hypothesis that (3.17) holds, concluding the proof of Lemma 3.26. ∎

According to Proposition 3.18, in order to complete the proof of Theorem 3.10, it remains only to analyze the case where the extreme point $s\in\{{0,1}\}^{m}$. Let $K=\{{k:s_{k}=1}\}$. Consider $W=\bigcap_{k\in K}\operatorname{ker}({\phi_{k}})$. Since $W$ is subcritical by hypothesis,

so $\operatorname{dim}({W})=0$, that is, $W=\{{0}\}$. Therefore the map $x\mapsto({\phi_{k}(x)})_{k\in K}$ from $V$ to the Cartesian product $\prod_{k\in K}V_{k}$ is injective.

since $s_{i}=0$ for all $i\notin K$. Thus it suffices to prove that

This is a special case of the following result.

Define $\Phi\colon V\to\prod_{k\in K}V_{k}$ by $\Phi(x)=({\phi_{k}(x)})_{k\in K}$. The hypothesis $\bigcap_{k\in K}\operatorname{ker}({\phi_{k}})=\{{0}\}$ is equivalent to $\Phi$ being injective. The product $\prod_{k\in K}\|f_{k}\|_{1}$ can be expanded as the sum of products

where the sum is taken over all $y=({y_{k}})_{k\in K}$ belonging to the Cartesian product $\prod_{k\in K}V_{k}$. The quantity of interest,

is likewise a sum of such products. Each term of the latter sum appears as a term of the former sum, and by virtue of the injectivity of $\Phi$, appears only once. Since all summands are nonnegative, the former sum is greater than or equal to the latter. Therefore

Having shown sufficiency for extreme points of $\mathcal{P}$, we apply Lemma 3.16 to conclude sufficiency for all $s\in\mathcal{P}$.

Communication lower bounds from Theorem 3.2

In addition to the concrete execution model and the geometric model, we use pseudocode in our examples. At a high level, these three different algorithm representations are related as follows:

A concrete execution is a sequence of instructions executed by the machine, according to the model detailed in Section 4.1. The concrete execution, unlike either the geometric model or pseudocode, contains explicit data movement operations between slow and fast memory.

The geometric model (Definition 2.2), is the most abstract, and is the foundation of the bounds in Section 3. Each instance (2.2) of the geometric model corresponds to a set of concrete executions, as detailed in Section 4.2.

The rest of this section is organized as follows. Section 4.1 describes the concrete execution model mentioned above, and Section 4.2 relates the concrete execution model to the geometric model. Section 4.3 states and proves the main communication lower bound result of this paper, Theorem 4.1. Section 4.4 presents a number of examples showing why the assumptions of Theorem 4.1 are in fact necessary to obtain a lower bound. Section 4.5 looks at one of these assumptions in more detail (“no loop splitting”), and shows that loop splitting can only improve (reduce) the lower bound. Finally, Section 4.6 discusses generalizations of the lower bound result to other machine models.

The hypothetical machine in our execution model has a two-level memory hierarchy: a slow memory of unbounded capacity and a fast memory that can store $M$ words (all data in our model have one-word width). Data movement between slow and fast memory is explicitly managed by software instructions (unlike a hardware cache), and data is copiedWe will use the word ‘copied’ but our analysis does not require that a copy remains, e.g., exclusive caches. at a one-word granularity (we will discuss spatial locality in Part 2 of this work). Every storage location in slow and fast memory (including the array elements $A_{j}(\phi_{j}(\mathcal{I}))$) has a unique memory address, called a variable; since the slow and fast memory address spaces (variable sets) are disjoint, we will distinguish between slow memory variables and fast memory variables. When a fast memory variable represents a copy of a slow memory variable $v$, we refer to $v$ as a cached slow memory variable; in this case, we assume we can always identify the corresponding fast memory variable given $v$, even if the copy is relocated to another fast memory location.

We define a sequential execution $E=({e_{1},e_{2},\ldots,e_{n}})$ as a sequence of statements $e_{i}$ of the following types:

$read(v)$: allocates a location in fast memory and copies variable $v$ from slow to fast memory.

$write(v)$: copies variable $v$ from fast to slow memory and deallocates the location in fast memory.

$compute(\{{v_{1},\ldots,v_{n}}\})$ is a statement accessing variables $v_{1},\ldots,v_{n}$.

$allocate(v)$ introduces variable $v$ in fast memory.

$free(v)$ removes variable $v$ from fast memory.

A sequential execution defines a total order on the statements, and thereby a natural notion of when one statement succeeds or precedes another. We say that a Read or Allocate statement and a subsequent Write or Free statement are paired if the same variable $v$ appears as an operand in both and there are no intervening Reads, Allocates, Writes, or Frees of $v$. A sequential execution is considered to be well formed if and only if

operands to Read (resp., Write) statements are uncached (resp., cached) slow memory variables,

operands to Allocate statements are uncached slow memory variables or fast memory variablesIf a fast memory variable $v$ in an $allocate(v)$ statement already stores a cached slow memory variable, then we assume the system will first relocate the cached variable to an available location within fast memory),

operands to Free statements are cached slow memory variables or fast memory variables,

every Read, Allocate, Write, and Free statement is paired, and

every Compute statement involving variable $v$ interposes between paired Read/Allocate and Write/Free statements of $v$, i.e., each operand $v$ in a Compute statements resides in fast memory before and after the statement.

Essentially, fast memory variables must be allocated and deallocated, either implicitly (Read/Write) or explicitly (Allocate/Free), while slow memory variables cannot be allocated/deallocated. Given the finite capacity of fast memory, we need an additional assumption to ensure the memory operations are well-defined. Given a well-formed sequential execution $E=({e_{1},\ldots,e_{n}})$, we define $footprint_{i}(E)$ to be the fast memory usage after executing statement $e_{i}$ in the program, i.e.,

Then $E$ is said to be $M$–fit (for fast memory size $M$) if $\max_{i}footprint_{i}(E)\leq M$.

It is of practical interest to permit variables to reside in fast memory before and after the execution, e.g., to handle parallel code as described in the next paragraph; however, this violates our notion of well-formedness. Rather than redefine well-formedness to account for this possibility, we take a simpler approach: given an execution that is well formed except for the $I$ (‘input’) and $O$ (‘output’) variables which reside in fast memory before and after the execution, we insert up to $M$ Reads and $M$ Writes at the beginning and end of the execution so that all memory statements are paired, and then later reduce the lower bound (on Reads/Writes) by $I+O\leq 2M$.

Although we will establish our bounds first for a sequential execution, they also apply to parallel executions as explained in Section 4.6. We define a parallel execution as a set of $P$ sequential executions, $\{{E_{1},E_{2},\ldots,E_{P}}\}$. In our parallel model, the global (‘slow’) memory for each processor is a subset of the union of the local (‘fast’) memories of the other processors. That is, for a given processor, each of its slow memory variables is really a fast memory variable for some other processor, and each of its fast memory variables is a slow memory variable for every other processor, unless it corresponds to a cached slow memory variable, in which case it is invisible to the other processors. (We could remove this last assumption by extending our model to distinguish between copies of a slow memory variable.) A parallel execution is well formed or (additionally) $\{{M_{1},\ldots,M_{P}}\}$–fit if each of its serial executions $E_{i}$ is well formed or (additionally) $M_{i}$–fit. Since well-formedness assumes no variables begin and end execution in fast/local memory, it seems impossible for there to be any nontrivial well-formed parallel execution. As mentioned above, we can always allow for a sequential execution with this property by inserting up to $I+O\leq 2M$ Reads/Writes, and later reducing the lower bound by this amount; we insert Reads/Writes in this manner to each sequential execution $E_{i}$ in a parallel execution.

2 Relation to the geometric model

Recall from Section 2 our geometric model (2.2):

Given an execution, we assume we can discern the expression $A_{j}(\phi_{j}(\mathcal{I}))$ from any variable which represents an array variable in the geometric model. The execution may contain additional variables that act as surrogates (or copies) of the variables specified in the program text. As an extreme example, an execution could use an array $A^{\prime}_{j}$ as a surrogate for the array $A_{j}$ in the computations, and then later set $A_{j}$ to $A^{\prime}_{j}$. In such examples, one can always associate each surrogate variable with the ‘master’ copy, and there is no loss of generality in our analysis to assume all variables are in fact the master copies.

We say a legal sequential execution of an instance of the geometric model is a sequential execution whose subsequence of Compute statements can be partitioned into contiguous chunks in one-to-one correspondence with $\mathcal{Z}$, and furthermore all $m$ array variables $A_{j}(\phi_{j}(\mathcal{I}))$ appear as operands in the chunk corresponding to $\mathcal{I}$. Given a possibly overlapping partition $\mathcal{Z}=\bigcup_{i=1}^{P}\mathcal{Z}_{i}$, a legal parallel execution is a parallel execution $\{{E_{1},\ldots,E_{p}}\}$ where each sequential execution $E_{i}$ is legal with respect to loop iterations $\mathcal{I}\in\mathcal{Z}_{i}$. Legality restricts the set of possible concrete executions we consider (for a given instance of the geometric model), and in general is a necessary requirement for the lower bound to hold for all concrete executions. For example, transforming the classical algorithm for matrix multiplication into Strassen’s algorithm (which can move asymptotically less data) is illegal, since it exploits the distributive property to interleave computations, and any resulting execution cannot be partitioned contiguously according to the original iteration space $\mathcal{Z}$. As another example, legality prevents loop splitting, an optimization which can invalidate the lower bound as discussed in Section 4.5.

We note that there are no assumptions about preserving dependencies in the original program or necessarily computing the correct answer. Restricting the set of executions further to ones that are correct may make the lower bound unattainable, but does not invalidate the bound.

3 Derivation of Communication Lower Bound

Now we present the more formal derivation of the communication lower bound, which was sketched in Section 2. The approach used here was introduced in [ITT04] and generalized in [BallardDemmelHoltzSchwartz11]. Here we generalize it again to deal with the more complicated algorithms considered in this paper.

Given an $M$–fit legal sequential execution $E=({e_{1},\ldots,e_{n}})$, we proceed as follows:

Break $E$ into $M$–Read/Write segments of consecutive statements, where each segment (except possibly the last one) contains exactly $M$ Reads and/or Writes. Each segment (except the last) ends with the $M^{\text{th}}$ Read/Write and the next segment starts with whatever statement follows. The last segment may have statements other than Reads/Writes at the end to complete the execution. (We will simply refer to these as Read/Write segments when $M$ is clear from context.)

Independently, break $E$ into Compute segments of consecutive statements so that the Compute statements within a segment correspond to the same iteration $\mathcal{I}$. (It will not matter that this does not uniquely define the first and last statements of a Compute segment.) Our assumption of a legal execution guarantees that there is one Compute segment per iteration $\mathcal{I}$. We associate each Compute segment with the (unique) Read/Write segment that contains the Compute segment’s first Compute statement.

Using the limited availability of data in any one Read/Write segment, we will use Theorem 3.2 to establish an upper bound $F$ on the number of complete Compute segments that can be executed during one Read/Write segment (see below). This is an upper bound on the number of complete loop iterations that can be executed.

Now, we can bound below the number of complete Read/Write segments by $\lfloor|\mathcal{Z}|/(F+1)\rfloor$. We add $1$ to $F$ to account for Compute segments that overlap two (or more) Read/Write segments. We need the floor function because the last Read/Write segment may not contain $M$ Reads/Writes. (Since we are doing asymptotic analysis, $|\mathcal{Z}|$ can often be replaced by the total number of Compute statements.)

Finally we bound below the total number of Reads/Writes by the lower bound on the number of complete Read/Write segments times the number of Reads/Writes per such segment, minus the number of Reads/Writes we inserted to account for variables residing in fast memory before/after the execution:

where we have applied our asymptotic assumption $F=\Omega(M)=\omega(1)$.

To determine an upper bound $F$, we will use Theorem 3.2 to bound the amount of (useful) computation that can be done given only $O(M)$ array variables. First, we discuss how to ensure that only a fixed number of array variables is available during a single Read/Write segment.

Given an $M$–fit legal sequential execution, consider any $M$–Read/Write segment. There are at most $M$ array variables in fast memory when the segment starts, at most $M$ array variables are read/written during the segment, and at most $M$ array variables remain in fast memory when the segment ends. If there are no Allocates of array variables, or if there are no Frees of array variables, then at most $2M$ distinct array variables appear in the segment (at most $M$ may already reside in fast memory, and at most $M$ more can be read or allocated). More generally, if there are no paired Allocate/Free statements of array variables, then at most $3M$ array variables appear in the segment (at most $M$ already reside in fast memory, at most $M$ can be read, and at most $M$ can be allocated). However, if we allow array variables to be allocated and subsequently freed, then it is possible to have an unbounded number of array variables contributing to computation in the same segment; this can occur in practice and we give concrete examples in the following section. Thus, we need an additional assumption in order to obtain a lower bound that is valid for all executions.

We will assume that the execution contains no paired Allocate/Free statements of array variables. However, we note that of these paired statements, we only need to avoid the ones where both statements occur within a given Read/Write segment; e.g., one could remove the preceding assumption by proving that at least $M$ Read/Write statements (of variables besides $v$) interpose every paired Allocate/Free of an array variable $v$. (This weaker assumption is equivalent to an assumption in [BallardDemmelHoltzSchwartz11, Section 2] that there are no ‘R2/D2 operands.’)

The communication lower bound is now a straightforward application of Theorem 3.2.

When $|\mathcal{Z}|=\Theta(F)$, i.e., the problem (iteration space) is not sufficiently large, the lower bound becomes $\Omega(M)-O(M)$, so the subtractive $O(M)$ term may dominate and lead to zero communication; this is increasingly likely as the ratio $|\mathcal{Z}|/F$ goes to zero. The parallel case also demonstrates this behavior in the ‘strong scaling’ limit, when the problem is decomposed to the point that each processor’s working set fits in their local memory (see Section 4.6). In the regime $|\mathcal{Z}|=O(F)$, a memory-independent lower bound [BDHLS12] provides more insight than the bound above. Let $W$ be the (unknown) number of Reads/Writes performed, and let $I$ and $O$ be the numbers of input/output variables residing in fast memory before/after execution. Then

4 Examples

We give four examples to demonstrate why our assumptions in Theorem 4.1 are necessary. Then, we discuss how one can sometimes deal with the presence of imperfectly nested loops, paired Allocate/Free statements, and an infeasible linear program to compute a useful lower bound; we give an example of this approach.

The following simple modification of matrix multiplication demonstrates how paired Allocate/Frees can invalidate our lower bound.

Suppose we know the arrays do not alias each other. Clearly $C$ only depends on the data when $i_{3}=N$, but we need to do all $N^{3}$ multiplications to compute $t$ correctly. But the same analysis from Section 3 applies to these loops as to matrix multiplication, suggesting a sequential communication lower bound of $\Omega(N^{3}/M^{1/2})$. However, it is clearly possible to execute all $N^{3}$ iterations moving only $O(N^{3}/M+N^{2})$ words, by hoisting the (unblocked) $i_{3}$ loop outside and blocking the $i_{1}$ and $i_{2}$ loops by $M/2$, doing $(M/2)^{2}$ multiplications in a Read/Write segment using $M/2$ entries each of $A(\cdot,i_{3})$ and $B(i_{3},\cdot)$, and (over)writing the values $C(i_{1},i_{2})$ to a single location in fast memory, which is repeatedly Allocated and Freed. Only when $i_{3}=N$ would $C(i_{1},i_{2})$ actually be written to slow memory. So in this case there are a total of $N^{3}-N^{2}$ paired Allocates/Frees, corresponding to the overwritten $C(i_{1},i_{2})$ operands. $\triangle$

This example also demonstrates how paired Allocate/Frees can invalidate our lower bound. Consider the following code:

Again, suppose we know that the arrays do not alias each other. Toward a lower bound, we ignore the initialization of $A$ (first loop nest) and only look at the second loop nest, a matrix multiplication. However, by computing entries of $A$ on-the-fly from $\mathcal{I}$ and discarding them (i.e., Allocating/Freeing them), one can beat the lower bound of $\Omega(N^{3}/M^{1/2})$ words for matrix multiplication. That is, by hoisting the (unblocked) $i_{2}$ loop outside and blocking the $i_{1}$ and $i_{3}$ loops by $M/2$, and finally writing each $A(i_{1},i_{3})$ to slow memory when $i_{2}=N$, we can instead move $O(N^{3}/M+N^{2})$ words. So, there are $N^{3}-N^{2}$ possible paired Allocate/Frees. $\triangle$

This example demonstrates how infeasibility of the linear constraints (3.1) of Theorem 3.2 can invalidate our lower bound. Consider the following code:

While infeasibility may be sufficient for there to be an unbounded number of array variables in a Read/Write segment, the previous two examples show that it is not necessary, since their linear programs are feasible. We will be more concrete about this relationship between infeasibility and unbounded data reuse in Part 2. $\triangle$

This example demonstrates how an execution which interleaves the inner loop bodies (an illegal execution) can invalidate our lower bound; see also Section 4.5. We will see in Section 4.5 that the lower bound for each split loop is no larger than the lower bound for the original loop. Consider splitting the two lines of the inner loop body in the Complicated Code example (see Section 1) into two disjoint loop nests (each over $\mathcal{I}\in\{{1,\ldots,N}\}^{6}$). We assume $\text{func}_{1}$ and $\text{func}_{2}$ do not modify their arguments, and that the arrays do not alias — the two lines share only read accesses to one array, $A_{3}$, so correctness is preserved. As will be seen later by using Theorem LABEL:thm6.1, the resulting two loop nests have lower bounds $\Omega(N^{6}/M^{3/2})$ and $\Omega(N^{6}/M^{2})$, resp., both better than the $\Omega(N^{6}/M^{8/7})$ of the original, and both these lower bounds are attainable. $\triangle$

Theorem 4.1 is enough for many direct linear algebra computations such as (dense or sparse) $LU$ decomposition, which do not have paired Allocate/Frees, but not all algorithms for the $QR$ decomposition or eigenvalue problems, which can potentially have large numbers of paired Allocates/Frees (see [BallardDemmelHoltzSchwartz11]). We can often deal with interleaving iterations, paired Allocates/Frees, and infeasibility of (3.1) by imposing Reads and Writes [BallardDemmelHoltzSchwartz11, Section 3.4]: we modify the algorithm to add (“impose”) Reads and Writes of array variables which are allocated/freed or repeatedly overwritten, apply the lower bound from Theorem 4.1, and then subtract the number of imposed Reads and Writes to get the final lower bound. (Note that we have already used a similar technique, above, to allow an execution to begin/end with a nonzero fast memory footprint.) We give an example of this approach (see also [BallardDemmelHoltzSchwartz11, Corollary 5.1]).

Consider computing $B=A^{k}$ using the following code, shown (for simplicity) for odd $k$ and initially $B=A$:

This example also illustrates that our results will only be of interest for sufficiently large problems, certainly where the floor function in the lower bound (4.1) is at least $1$. $\triangle$

The above approach covers many but not all algorithms of interest. We refer to reader to [BallardDemmelHoltzSchwartz11, Sections 3.4 and 5] for more examples of imposing Reads and Writes, and [BallardDemmelHoltzSchwartz11, Section 4] on orthogonal matrix factorizations for an important class of algorithms where a subtler analysis is required to deal with paired Allocates/Frees.

Imposing Reads and Writes may fundamentally alter the program, so the lower bound obtained for the modified code need not apply to the original code. In the example above, one could reorderFor simplicity, and without reducing data movement, this code performs $(k-1)N^{2}$ additional $+$ operations. the original code to

5 Loop Splitting Can Only Help

Here we show that loop splitting can only reduce (improve) the communication lower bound expressed in Theorem 4.1. More formally, we state this as the following.

While loop splitting appears to always be worth attempting, in practice data dependencies limit our ability to perform this optimization; we will discuss the practical aspects of loop splitting further in Part 2 of this work.

6 Generalizing the machine model

Earlier we said the reader could think of a sequential algorithm where the fast memory consists of a cache of $M$ words, and slow memory is the main memory. In fact, the result can be extended to the following situations:

If a sequential machine has a memory hierarchy, i.e., multiple levels of cache (most do), where data may only move between adjacent levels, and arithmetic done only on the “top” level, then it is of interest to bound the data transferred between every pair of adjacent levels, say $i$ and $i+1$, where $i$ is higher (faster and closer to the arithmetic unit) than $i+1$. In this case we apply our model with $M$ representing the total memory available in levels 1 through $i$, typically an increasing function of $i$.

One may also ask what value of $M$ to use for each processor. Suppose that each processor has $M_{\text{max}}$ words of fast memory, and that the total problem size of all the array entries accessed is $M_{\text{array}}$. So if each processor gets an equal share of the data we use $M=M_{\text{array}}/P\leq M_{\text{max}}$. But the lower bound may still apply, and be smaller, if $M$ is larger than $M_{\text{array}}/P$ (but at most $M_{\text{max}}$). In some cases algorithms are known that attain these smaller lower bounds (e.g., matrix multiplication in [2.5D_EuroPar]), i.e., replicating data can reduce communication.

In Section 7.3, we discuss attainability of these parallel lower bounds, and reducing communication by replicating data.

The simplest possible hierarchical machine is the sequential one with multiple levels of memory discussed above. But real parallel machines are similar: each processor has its own memory organized in a hierarchy. So just as we applied our lower bound to measure memory traffic between levels $i$ and $i+1$ of cache on a sequential machine, we can similarly analyze the memory hierarchy on each processor in a parallel machine.

Finally, people are building heterogeneous parallel machines, where the various processors, memories, and interconnects can have different speeds or sizes. Since minimizing the total running time means minimizing the time when the last processor finishes, it may no longer make sense to assign an equal fraction $|\mathcal{Z}|/P$ of the work and equal subset of memory $M$ to each processor. Since our lower bounds apply to each processor independently, they can be used to formulate an optimization problem that will give the optimal amount of work to assign to each processor [Hetero_SPAA11].

(Un)decidability of the communication lower bound

In Section 3, we proved Theorem 3.2, which tells us that the exponents $s\in^{m}$ satisfy the inequalities (3.1), i.e.,

The proof of Theorem 5.1 is built upon several smaller results.

Generate a list of all nonempty subsets of $V$ having at most $\operatorname{dim}({V})$ elements. Test each subset for linear independence, and discard all which fail to be independent. Output a list of those which remain. ∎

We do not require this list to be free of redundancies.

To the given family of inequalities, adjoin the $2m$ inequalities $s_{j}\geq 0$ and $-s_{j}\geq-1$. $\mathcal{P}$ is the convex polytope defined by all inequalities in the resulting enlarged family. Express these inequalities as $\langle s,v_{\alpha}\rangle\geq c_{\alpha}$ for all $\alpha\in A$, where $A$ is a finite nonempty index set.

Create a list of all subsets $B\subset A$ with cardinality equal to $m$. There are finitely many such sets, since $A$ itself is finite. Delete each one for which $\{{v_{\beta}:\beta\in B}\}$ is not linearly independent. For each subset $B$ not deleted, compute the unique solution $\tau$ of the system of equations $\langle\tau,v_{\beta}\rangle=c_{\beta}$ for all $\beta\in B$. Include $\tau$ in the list of all extreme points, if and only if $\tau$ satisfies $\langle\tau,v_{\alpha}\rangle\geq c_{\alpha}$ for all $\alpha\in A\setminus B$. ∎

There exists an algorithm which takes as input a vector space HBL datum $({V,({V_{j}}),({\phi_{j}})})$, an element $t\in^{m}$, and a subspace $\{{0}\}<W<V$ which is critical with respect to $t$, and determines whether $t\in\mathcal{P}({V,({V_{j}}),({\phi_{j}})})$.

Theorem 5.1 and Proposition 5.5 will be proved inductively in tandem, according to the following induction scheme. The proof of Theorem 5.1 for HBL data in which $V$ has dimension $n$, will rely on Proposition 5.5 for HBL data in which $V$ has dimension $n$. The proof of Proposition 5.5 for HBL data in which $V$ has dimension $n$ and there are $m$ subspaces $V_{j}$, will rely on Proposition 5.5 for HBL data in which $V$ has dimension strictly less than $n$, on Theorem 5.1 for HBL data in which $V$ has dimension strictly less than $n$, and also on Theorem 5.1 for HBL data in which $V$ has dimension $n$ and the number of subspaces $V_{j}$ is strictly less than $m$. Thus there is no circularity in the reasoning.

Let $({V,({V_{j}}),({\phi_{j}})})$ and $t,W$ be given. Following Notation 3.23, consider the two HBL data ({W,({\phi_{j}(W)}),({\phi_{j}\big{|}_{W}})}) and $({V/W,({V_{j}/\phi_{j}(W)}),({[\phi_{j}]})})$, where $[\phi_{j}]\colon V/W\to V_{j}/\phi_{j}(W)$ are the quotient maps. From a basis for $V$, bases for $V_{j}$, a basis for $W$, and corresponding matrix representations of $\phi_{j}$, it is possible to compute the dimensions of, and bases for, $V/W$ and $V_{j}/\phi_{j}(W)$, via row operations on matrices. According to Lemma 3.25, $t\in\mathcal{P}({V,({V_{j}}),({\phi_{j}})})$ if and only if t\in\mathcal{P}({W,({\phi_{j}(W)}),({\phi_{j}\big{|}_{W}})})\cap\mathcal{P}({V/W,({V_{j}/\phi_{j}(W)}),({[\phi_{j}]})}).